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Practical application of flow boiling to ground- and space-based thermal management 
systems hinges on the ability to predict the system’s heat removal capabilities under 
expected operating conditions. Research in this field has shown that the heat transfer 
coefficient within two-phase heat exchangers can be largely dependent on the experienced 
flow regime. This finding has inspired an effort to develop mechanistic heat transfer 
models for each flow pattern which are likely to outperform traditional empirical 
correlations. As a contribution to the effort, this work aimed to identify the heat transfer 
mechanisms for the slug flow regime through analysis of individual Taylor bubbles. 

An experimental apparatus was developed to inject single vapor Taylor bubbles into co- 
currently flowing liquid HFE 7100. The heat transfer was measured as the bubble rose 
through a 6 mm inner diameter heated tube using an infrared thermography technique. 
High-speed flow visualization was obtained and the bubble film thickness measured in an 
adiabatic section. Experiments were conducted at various liquid mass fluxes (43-200 


kg/m’s) and gravity levels (0.01g-1.8g) to characterize the effect of bubble drift velocity 


on the heat transfer mechanisms. Variable gravity testing was conducted during a NASA 
parabolic flight campaign. 

Results from the experiments showed that the drift velocity strongly affects the 
hydrodynamics and heat transfer of single elongated bubbles. At low gravity levels, 
bubbles exhibited shapes characteristic of capillary flows and the heat transfer 
enhancement due to the bubble was dominated by conduction through the thin film. At 
moderate to high gravity, traditional Taylor bubbles provided small values of enhancement 
within the film, but large peaks in the wake heat transfer occurred due to turbulent vortices 
induced by the film plunging into the trailing liquid slug. Characteristics of the wake heat 
transfer profiles were analyzed and related to the predicted velocity field. Results were 
compared and shown to agree with numerical simulations of colleagues from EPFL, 
Switzerland. 

In addition, a preliminary study was completed on the effect of a Taylor bubble passing 
through nucleate flow boiling, showing that the thinning thermal boundary layer within the 


film suppressed nucleation, thereby decreasing the heat transfer coefficient. 
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CHAPTER 1: INTRODUCTION 


1.1 BACKGROUND/MOTIVATION 


Flow boiling is a type of heat transfer that combines convective and phase-change 
characteristics of the working liquid to improve the efficiency of thermal management 
systems. Traditional single-phase heat exchangers often require high liquid flow rates, and 
subsequently large pumping power, to effectively cool high power density loads. The 
addition of boiling to the flow provides significant heat transfer enhancement due to the 
latent heat of the fluid and convective forces contributed by the vapor. Depending on the 
mass flow rate, vapor quality (defined as the mass fraction of vapor to total mass), and wall 
heat flux applied, the flow may assume several regimes. The most common flow patterns 
for vertical-upward flow boiling in a tube are illustrated in Figure 1.1. Bubbly flow (Figure 
1.1a) is characterized by small bubbles which are nucleated from the heated tube wall and 
travel upward with the flow. As more vapor (bubbles) are added downstream, the small 
bubbles coalesce to form larger vapor plugs that are nearly the diameter of the tube and are 
separated by liquid slugs. This regime is classified as slug flow (Figure 1.1b). Further 
increase of the vapor quality leads to churn flow (Figure 1.1c) in which the individual plugs 
break down into a chaotic and oscillatory column of vapor surrounded by liquid on the 
wall. Annular flow (Figure 1.1d) is obtained when the shear force exerted by the vapor core 


becomes strong enough to maintain a steady liquid film on the tube wall. 


a/g 


Flow 


a) b) c) d) 


Figure 1.1: Schematic of common flow boiling regimes for vertical-upward flow within a heated 
tube: a) bubbly flow, b) slug flow, c) churn flow, and d) annular flow. 


The study of flow boiling under reduced gravity conditions has gained interest in recent 
decades due to the growing heat dissipation requirements placed on space-based heat 
exchanger systems. As one of its top technical challenges, NASA has targeted two-phase 
heat transfer loops as a way to more efficiently transfer large amounts of heat with little 
temperature change of the thermal load [1]. An effort to advance the technology readiness 
level of these systems to TRL 6 will require a more thorough understanding of heat transfer 
mechanisms in partial or microgravity (a/g < 0.01) environments. Thermal control 
systems can then be designed to utilize two-phase loops as a means to dissipate heat from 
high power density electronics, while reducing launch weight and volume. To acquire the 
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necessary understanding and expertise in the field, NASA and other national space 
agencies have looked to two-phase researchers to study the effect of gravity on flow boiling 


heat transfer. 


1.2 MICROGRAVITY RESEARCH 


One of the earliest studies on flow boiling under reduced gravity conditions was by Saito 
et al. [2]. Their work analyzed the flow of water over a heated rod oriented parallel to the 
length of a square channel and horizontally with respect to the ground (or aircraft floor). 
Data collected during parabolic flights illustrated larger bubble departure diameters during 
the microgravity phase due to the reduction of buoyancy force acting on the bubbles. This 
acted to decrease the heat transfer coefficient on the underside of the heated rod and 
increase it on the top. The variation in heat transfer coefficient, however, was fairly small 
in comparison to the large differences in flow regime. 

Around the same time, Lui et al. [3] carried out flow boiling heat transfer experiments 
with subcooled R113 in a 12 mm ID tubular test section. Heat transfer coefficients were 
found to increase by 5 to 20% under microgravity conditions, typically increasing with 
quality. Turbulence caused by greater movement of vapor was given as the reason for the 
trend. Contradictory results to Lui et al. were found by both Fore et al. [4] and Rite and 
Rezkallah [5] who conducted experiments on bubbly and slug flows in tubular channels. A 
decrease in heat transfer coefficient by as much as 50% at microgravity was seen by Rite 
and Rezkallah at low gas qualities, but the discrepancy between microgravity and normal 
gravity heat transfer became much smaller as the quality neared annular flow transition. 
Fore et al. also found lower heat transfer coefficients in microgravity slug and bubbly flow. 


In both cases, the decrease in heat transfer at low qualities under microgravity conditions 
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was attributed to the reduction in turbulence behind bubbles, resulting from a smaller 
relative velocity between the bubbles and liquid (the drift velocity) in the absence of 
buoyancy force. 

Ohta [6] also found the effect of flow dynamics important in his study of R113 
evaporation in an 8 mm ID tubular channel. Unlike Fore et al. and Rite and Rezkallah, Ohta 
observed little difference in heat transfer coefficient for bubbly flow despite much larger 
microgravity bubble departure diameters at low flow rates (G=150 kg/m/’s), as can be seen 
in Figure 1.2. He asserted that this confirmed the importance of turbulence near the heated 
wall generated by bubble nucleation, rather than the bubble behavior in the central flow. 
When nucleation was suppressed, however, as in moderate vapor quality and low to 
moderate heat fluxes, a variation in heat transfer coefficient with gravity was observed: it 
was enhanced in high gravity (2g) and decreased in microgravity. It was explained that the 
frequency and length of the passing disturbance waves in the annular flow pattern increased 
in the 2g period and decreased in microgravity, leading to a thinning and thickening of the 
liquid film thickness, respectively. When the wall heat flux was increased under the same 
quality conditions to the extent that nucleation occurred within the liquid film, Ohta found 


that the heat transfer coefficient was again unaffected by gravity level. 
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Figure 1.2: Time trace of heat transfer coefficient and acceleration with corresponding flow 
visualization images for G=150 kg/m’s. Figure from Ohta [6] (used with permission). 


The group of Celata and Zummo have completed extensive work on the characterization 
of microgravity heat transfer and the determination of gravity’s mechanistic effects. 
Utilizing a Pyrex tube heated using a helical resistor, visualization and heat transfer 
measurements were obtained for boiling FC-72. Early results [7] found that for subcooled 
boiling at low mass fluxes, a decrease in heat transfer was experienced at microgravity, 
while for high mass fluxes, no variation was seen. Additionally, it was observed that 
regardless of mass flux, vapor qualities higher than 0.3 saw little change in heat transfer 
coefficient between normal and microgravity. The authors suggest that a growth in bubble 
size in microgravity is typically accompanied by a difference in heat transfer, while similar 
bubble sizes and distribution lead to negligible difference. 

More recently, Zummo et al. [8] found that the variation in tube wall temperature and 
heat transfer coefficient during the microgravity period depended on the axial measurement 


location in the tube. They also reiterated their previous findings that above a critical value 


of fluid velocity, inertial forces dominate over surface tension and buoyancy forces acting 
on the bubbles, leading to little change in bubble size with gravity variation. As a result, 
the turbulence created by the bubbles that disrupt the thermal boundary layer will be similar 
for all gravity levels, and one should expect that the wall temperature will not be affected. 

In agreement with many of the previous works, Narcy et al. [9] found that at high mass 
fluxes (G > 400 kg/m’s) the flow regime and heat transfer coefficients were unchanged by 
a variation in gravity level. When the mass flux was reduced to G = 200 kg/m’s, the heat 
transfer coefficient seemed to be smaller under microgravity conditions than in normal 
gravity, but further work was suggested. Their experiment consisted of a6 mm ID sapphire 
tube heated resistively with a thin indium tin oxide (ITO) coating which allowed for 
simultaneous heat transfer measurements and high speed flow visualization. After 
collecting more data, Narcy et al. [10] concluded that for mass fluxes between G = 100 
kg/m?s and G = 400 kg/m’s, saturated boiling heat transfer is weakly affected by gravity 
level. However, they found that at low heat flux and subcooled conditions, the heat transfer 
coefficient is 20% lower while in a microgravity environment. They theorized that the 
degradation is caused by lower bubble formation frequency in microgravity. 

Narcy and Colin [11] provided further comparison between normal and microgravity 
heat transfer, noting that at G = 200 kg/m?s and qualities smaller than 0.15, where bubbly 
and slug flow regimes were present, microgravity heat transfer coefficients were smaller 
than those in normal gravity. Regardless of quality and flow regime, it was found that a 
mass flux of G = 50 kg/m?s always yielded a lower heat transfer coefficient in microgravity. 
They suggested that the nucleate boiling regime remained dominant even when annular 


flow was present. 


While not explicitly commenting the effect of gravity on heat transfer, Westheimer and 
Peterson [12] found that smaller heat fluxes were needed to cause flow regime transitions 
in reduced gravity environments. Through observations of boiling R113 in a glass annular 
heat exchanger, they showed that a flow regime characteristic of higher vapor quality in 
normal gravity could be seen at a lower quality in microgravity, while maintaining constant 
flow and heating conditions. This trend has been partly noted in work described above and 
plays an important role in this proposed study. Specifically, it has been shown in recent 
years that the variation of heat transfer coefficient with gravity level is likely attributed to 
differences in flow dynamics and flow regime. This, then, highlights the necessity for two 
main areas of research important for the prediction of gravity’s effect: reduced gravity flow 
regime maps and mechanistic heat transfer studies for each flow regime. Numerous studies 
have been conducted to create accurate flow pattern maps under different gravity 
conditions, and results are strongly influenced by the fluid used as well as the experimental 
setup. However, of greater need is the determination of the specific heat transfer 
mechanisms at work for the flow regimes experienced in microgravity conditions. Slug 


flow is of particular interest for this study. 


1.3 SINGLE ELONGATED BUBBLE RESEARCH 


Slug flow is characterized by elongated bubbles separated by liquid slugs and can be 
broken down into “unit cells”, each consisting of an individual bubble as shown in Figure 
1.3. For the vertical, upward co-current flow configuration, the bubble rises with a velocity 
(U;) greater than the liquid (U,) due to buoyancy and the non-uniform velocity profile 
present in the tube. The relative velocity between the bubble and liquid can be defined as 


the drift velocity, Ug = U, — U;. In order to replace the volume left behind the bubble as 
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it rises, liquid from the leading slug falls around the bubble, forming a film with thickness 
0. The film falls due to gravity and is balanced by the shear force experienced at the tube 
wall. At the tail, the film plunges as a circular wall jet into the trailing slug creating vortices 
whose complexity depend on Uq. The bubble shape is dictated by the ratio of buoyancy to 


. F - D2 
surface tension forces as described by the Bond number, Bo = os 2 
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Figure 1.3: Diagram describing the liquid flow pattern around a rising Taylor bubble. 


In detailed analysis of the flow field around a bubble, two reference frames can be used: 
a laboratory reference frame where the bubble rises in a stationary tube and a bubble 
reference frame which moves with the bubble as shown in Figure 1.4. This assumes the 
bubble does not grow as it rises in the tube. In this dissertation, unprimed and primed 
quantities refer to area average velocities in the laboratory and bubble reference frames, 


respectively. It should be noted that the drift velocity, Ug = U, — U,; = U; — Uj, is 
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independent of the chosen reference frame as is the plunging velocity, U, = U; + U; = 


U; — Uj. The average film velocity, Uy, can be calculated from a mass balance and the film 


thickness: 
(U, — Up) R? 
U,; = ————- - 1 
ff (2R5—62) ~? o 
Laboratory Reference Frame Bubble Reference Frame 


Figure 1.4: Diagram of velocity profiles ahead of the bubble and in the near wake region in the 
laboratory (left) and bubble (right) reference frames. 


1.3.1 Taylor Bubble Hydrodynamics 

Bubbles with Bo > 40 are considered macroscale [13], exhibit bullet-shaped profiles, and 
are commonly described as Taylor bubbles (from the early work by Davies and Taylor 
[14]). Bubbles of this type have been widely studied; a review of the literature prior to 1992 
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is given in Fabre and Line [15], which summarizes fundamental aspects of the bubble 
hydrodynamics. Nicklin et al. [16] proposed a correlation for bubble velocity, Up, of the 
form, 

Up = CU, + Up (2) 
where C is a constant that defines the bubble rise velocity depending on the fluid velocity 
profile, and Up,o is the bubble rise velocity in a stagnant fluid column. Values of C were 
theoretically and experimentally determined to be 2 and 1.2 for a Taylor bubble flowing 
within fully developed laminar and turbulent flows, respectively [16,17]. White and 
Beardmore [18] proposed a graphical relation for Up,o based on experimental measurements 
of air bubbles in a variety of stagnant liquids. Results were broken into several regimes 
depending on the relative influence of inertial, surface tension, and viscous forces acting 
on the bubble. Brown [19] modified the potential flow solution for the flow field around 
rising Taylor bubbles to account for liquid viscosity and developed a correlation that agreed 
well with experimental results when viscous forces could be ignored at the bubble nose. 
Viana et al. [13] combined bubble velocity data from the literature with their own data to 
create a universal correlation for Up,o based on Bo and a buoyancy Reynolds number similar 
to the Grashof number, Gr. Rattner and Garimella [20] recently completed an experimental 
study on intermediate scale slug flows and proposed an extension of Bendiksen’s [21] 
prediction for Up,o. Polonsky et al. [22] analyzed high speed video of rising Taylor bubbles 
in air-water vertical flow to characterize the bubble velocity and shape. Oscillations in the 
bubble tail were found to be of higher frequency and amplitude for longer bubbles. 

As the bubble rises, a liquid film is formed between the liquid/vapor interface and the 


tube wall that thins as it moves towards the tail, eventually reaching a constant thickness 
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(6) where the liquid velocity profile within the film has become fully developed. Campos 
and Carvalho [23] proposed an equation (eqn. (3)) to calculate the distance, Z, at which the 


film stabilizes. 
a (evo) (3) 


Nogueira et al. [24] compared experimentally measured values of Z to predictions by 


Campos and Carvalho and found that the agreement was heavily dependent on the 
measured film Reynolds number, Rey Pi “er where Uf is the mean velocity in the fully 
developed film. Z was underpredicted at film Reynolds numbers below 40, but 
overpredicted above Rey 7 = 80. 

Prediction of the fully developed film thickness has been the subject of both theoretical 
and experimental work. Goldsmith and Mason [25] developed a relation for 6 by solving 


the Navier-Stokes equations for a laminar thin film around a Taylor bubble. The thin film 


assumption was later relaxed by Brown [19] to obtain a similar result (eqn. (4)). 


v / 
b= Ls ((R — 6)2U, — RU) ; (4) 


Llewellin et al. [26] approximated 6Ofor air bubbles rising in a variety of Newtonian 
fluids by correlating the bubble volume to its length. It was found that the relative film 
thickness (6/R) was solely a function of the non-dimensional parameter Ny = gp" (this is 
also equal to the square root of the Galileo number, Ga). An empirical correlation was 


proposed based on their collected data as well as the data of Nogueira et al. [24], 


6 
Sat b tanh(c —dlogio Nr) (5) 


11 


where the constant values are a=0.204, b=0.123, c=2.66, and d=1.15. Rattner and 
Garimella [20] measured film thickness, bubble velocity, and void fraction for Taylor 
bubble trains in 6, 8, and 9.5 mm ID tubes using an optical flow visualization facility. 
Bubble film thicknesses were found to be relatively independent of the liquid flow rates, 
and depended primarily on the tube diameter. A series of equations was proposed to predict 
the film thickness for intermediate scale Taylor flows. 

At the tail, the liquid film plunges into the liquid slug trailing the bubble creating 
vortices in the wake whose structures vary in complexity from organized toroids to chaotic 
shedding. Characterization of the flow field in this region has seen significant progress in 
recent years due in large part to the advancement of experimental techniques. Campos and 
Carvalho [23] performed a photographic study of rising air bubbles in stagnant 
water/glycerol mixtures seeded with dye to determine the structure of the bubble wake. 
The parameter Ny was varied by adjusting the ratio of water to glycerol. Three wake 
patterns were identified as shown in Figure 1.5. For Nr < 500 (Figure 1.5a) the wake 
exhibited a laminar, closed, axisymmetric, toroidal vortex with internal recirculatory flow. 
A transition region occurred for 500 < N¢ < 1500 (Figure 1.5b) where the toroidal vortex 
remained intact, but became non-axisymmetric due to oscillations in the bubble tail. For 
Nr > 1500 (Figure 1.5c), the vortices in the wake were shed randomly, exhibiting 


turbulent characteristics that decayed as they moved downstream. 
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a) b) Cc) 


Figure 1.5: Wake patterns described by Campos and Carvalho [23]: a) laminar, closed, 
axisymmetric, toroidal vortex with internal recirculatory flow; b) asymmetric, closed, toroidal 
vortex; c) open, turbulent vortices. Used with permission. 


Liu et al. [27] compared these boundaries to wake velocity profiles for vapor bubbles 
rising in stagnant liquid nitrogen. Results from their particle image velocimetry (PIV) 
technique revealed that due the significant property differences between water (used by 
Campos and Carvalho) and liquid nitrogen, the prediction of wake regimes was not 
applicable to cryogenics. Transitional and laminar wake patterns were observed for Ny > 


1500. Pinto et al. [28] redefined Campos and Carvalho’s criteria in terms of a Reynolds 


number based on the relative velocity between the bubble and liquid, Rey, = p eae They 


proposed that the wake would be fully turbulent for Rey, > 525 in the case of liquid 
flowing co-currently with the bubble, but noted that the limit had not been experimentally 


confirmed. 


13 


Kawaji et al. [29] observed turbulent wake characteristics and suggested the vortices 
were due to a Helmholtz-type instability resulting from the relative velocity difference 
between the plunging film and the trailing slug. They found that the distance behind the 
bubble at which the film became unstable (the “penetration length”) remained relatively 
constant with respect to bubble length once the flow within the film became fully 
developed. 

Shemer et al. [30] used PIV to study the velocity profiles in the wake of single bubbles 
rising in both laminar and turbulent background flows. They found the wake flow field 
could be effectively turbulent even when the base flow was laminar. Turbulence quantities 
calculated from the data showed that the initial mixing process occurred in the near wake 


region and persisted a few diameters downstream of the bubble tail. 


1.3.2 Taylor Bubble Heat Transfer 

Hetsroni and Rozenblit [31] briefly touched upon the mechanisms of heat transfer in slug 
flow. A thin heated film was installed within a 74 mm ID tube and coated with black paint 
to allow temperature measurement using an IR camera. The wall temperature during the 
passage of a Taylor bubble was observed to be higher compared to the liquid slug behind 
the bubble. This suggests a higher heat transfer coefficient in the liquid slug than in the 
liquid film, but no quantitative data was presented. 

Using a similar infrared technique, Babin et al. [32-34] investigated the wake heat 
transfer enhancement for single air bubbles rising in stagnant and vertically flowing water 
in 26 and 44 mm ID tubes. They observed that the heat transfer coefficient rose rapidly just 
behind the bubble tail before decaying to the original single-phase value several hundred 


diameters downstream. It was also seen that bubbles moving co-currently with turbulent 


14 


background flow created smaller heat transfer enhancement compared to those moving in 
laminar background flow. Heat transfer measurements were compared to PIV velocity 
measurements obtained by Shemer et al. [30] to show that the growth and decay of velocity 
fluctuations in the wake may correspond to the heat transfer coefficient profile observed 
for turbulent background flows. 

A model for macroscale gas/liquid slug flow was developed analytically by Barnea and 
Yacoub [35] in which an energy balance was performed on fluid cross-sections to 
determine the local wall and fluid temperatures as a function of time. For the case of 
constant heat flux, the wall temperature was found to be higher during the passing of a 
liquid slug compared to the liquid film. While this suggests higher heat transfer in the film, 
it should be noted that the model is based on the absorptive and convective capabilities of 


each cross-section and not on observed mechanisms. 


1.3.3 Capillary Bubble Heat Transfer 

The advent of microscale (defined as Bo < 0.9 — 19.7 by [18,21,36—38]) heat exchangers 
has spurred work into the behavior of microchannel slug flow. Slug flow under these 
conditions is characterized by elongated bubbles, also known as capillary bubbles, with 
near-spherical caps at the nose and the tail. Several modeling efforts and numerical studies 
of local heat transfer around capillary bubbles have been completed. A model by Jacobi 
and Thome [39] featured a two-zone representation of evaporation for elongated bubble 
flows where the regime was divided into a liquid slug region and a thin film region. It was 
suggested that the main mechanism of heat transfer was evaporation of the thin film trapped 
between the bubble and heated channel wall. This model was modified by Thome et al. 


[40] to include three zones as shown in Figure 1.6: a liquid slug, an evaporating elongated 
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bubble, and a vapor slug. As with the two-zone model, evaporation in the thin liquid film 
provided heat transfer several orders of magnitude higher than the single phase heat transfer 


due to the liquid slug. 
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Figure 1.6: Diagram of the three-zone heat transfer model for the elongated bubble regime in 
microchannels, illustrating a triplet comprised of a liquid slug, an elongated bubble, and a vapor 
slug (from Thome et al. [40], used with permission). 


Magnini et al. [41] numerically obtained the shape, length, and local heat transfer of 
Taylor bubbles during boiling of several fluids in a 0.5 mm circular channel. They found 
that as the bubble entered a heated channel containing a developing thermal boundary 
layer, evaporation of the liquid film removed heat from the fluid and caused the heat 
transfer to become larger than for single-phase flow. The heat transfer coefficient rose 
monotonically from the bubble nose towards the tail, with the highest values occurring in 
the bubble wake region. Based on these results, they modified the three-zone model by 
Thome et al. [40] to include unsteady conduction through the liquid, and obtained better 


agreement with the simulations. 
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Extending their previous study, Magnini et al. [42] analyzed the effect of leading and 
sequential bubbles on evaporating slug flow heat transfer of refrigerants in a 0.5 mm tube. 
Two bubbles were injected into a thermally-developing flow with constant heat flux 
boundary conditions applied to the tube wall. It was seen that the growth of the trailing 
bubble accelerated the leading bubble, but not the liquid within the film around the leading 
bubble. This increased the difference between liquid film velocity and bubble vapor 
velocity, thereby promoting local instability at the interface. It is also observed that the 
average heat transfer for the trailing bubble is nearly twice that observed for single-phase 
flow and 60% higher than the leading bubble, a result attributed to the re-development of 
the thermal boundary layer after the first bubble transit and perturbations of the second 
bubble. The authors created a new model for the heat transfer coefficient, breaking the flow 
structure into an adherent film, present in the bubble region and the liquid slug, and a 
recirculating flow within the slug as illustrated in Figure 1.7. A comparison to numerical 
simulations yielded good agreement when superheating of the liquid slug is taken into 
account. Many other numerical studies have been conducted on both the hydrodynamics 
and heat transfer of slug flow in microchannels and have been summarized in a review by 


Talimi et al. [43]. 
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Figure 1.7: Schematic of microchannel slug flow decomposition into adherent film and 
recirculation region within the liquid slug (from Magnini et al. [42], used with permission). 


The experimental heat transfer work to date has largely focused on measuring the 
overall heat transfer enhancement of slug trains with respect to single phase flow, rather 
than on the mechanisms of heat transfer around each bubble. For example, Walsh et al. 
[44,45] utilized an IR technique in which the outer wall temperature of a 1.5 mm ID 
stainless steel tube was measured when air-water bubble trains were present. The observed 
wall temperature was used as a boundary condition for a thermal resistance problem to 
obtain the time-averaged heat transfer coefficient. The maximum heat transfer 
enhancement over fully developed Poiseuille flow occurred at a liquid slug length to 
diameter ratio of unity. A correlation to predict the fully developed Nusselt number at other 
ratios was proposed. 

Mehta and Khandekar [46] utilized a similar IR technique in bubble train experiments 
using deionized water and air within a5 mm x 5 mm square minichannel. The heat transfer 
coefficient along the channel length was calculated using the experimentally measured 
channel wall temperature contours. An enhancement in heat transfer coefficient of 1.2 to 2 
times over thermally developing single-phase flow was observed depending on the axial 


location. 
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1.4 OBJECTIVES 


The objective of this work was to identify the contributions of various heat transfer 
mechanisms acting on a heated tube in the presence of a single elongated bubble at a variety 
of gravitational accelerations. Additionally, observations of the bubble shape and dynamics 
were made to allow for predictions of the flow field, which complemented the heat transfer 
measurements. The collected data can be utilized to guide mechanistic modeling of slug 
flow heat transfer and to provide validation for their results. Successful models may be 
used as design tools for future earth- and space-based two-phase heat exchanger systems, 


thereby satisfying NASA’s goal of this research. 
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CHAPTER 2: EXPERIMENT DESCRIPTION 


To characterize the heat transfer and dynamics of rising Taylor bubbles, a flow boiling 
experiment was conducted in which measurements of the local wall heat transfer and film 
thickness along with high speed images were obtained as single bubbles of varying length 
rose in a vertical column containing upward flowing liquid. The study of single bubbles 
was chosen in lieu of bubble trains so flow conditions upstream and downstream of the 
bubbles could be measured, reducing the complexity in approximating the flow patterns 
and understanding the heat transfer profile around each bubble. 

A flow boiling test rig developed in the Phase Change Heat Transfer Lab was utilized 
for this study. The rig, which was designed for parabolic flight testing, includes the flow 
loop and various electronics required for testing. Figure 2.1 shows the rig in flight 
configuration with a secondary containment chamber surrounding the flow loop to avoid 
potential fluid leakage into the aircraft cabin. Three computers and four monitors were used 
to control the flow loop and acquire infrared, high-speed visual, and film thickness data. 
Several layers of safety features (overheat protection, circuit protection, strength 
requirements) were included per NASA specifications. The flow loop and infrared camera 


technique are described in the following sections. 
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Figure 2.1: Experimental rig in parabolic flight configuration. 


2.1 FLOW LOOP 


A schematic of the experiment flow loop is shown in Figure 2.2. The working fluid was 
3M Novec HFE 7100 (C4F90CHs), a non-toxic, dielectric fluid with a normal boiling 
temperature of 60°C. Properties at saturation conditions for 1 bar of pressure are 


summarized in Table 2.1. 


Table 2.1: Summary of HFE 7100 properties at saturations conditions for 1 bar. 


Property Saturation Value 
Tsat PC) 60 

p, [kg/m?] 1372 

Ay [MJ/kg] 112.1 

Lt, [cP] 0.375 

o, [N/m] 0.128 
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Figure 2.2: Schematic of rig flow loop. 


HFE 7100 was pumped in a subcooled liquid state using a gear pump (Micropump 


L21755) as the flow rate was measured by a turbine flowmeter (Omega FLR1009). The 


liquid was heated to near saturation at the test section inlet using a stainless steel preheater 


powered by a modified 1000W computer power unit (Silverstone SST-ST1000-P) and 


controlled using pulse width modulation via a LabVIEW™ interface. 


The fluid then entered a section of the flow loop designed to create and release single 


Taylor bubbles. This section consisted of a bubble generation segment and a bypass 


segment connected at the downstream end by a three-way make-before-break valve as 


shown in Figure 2.3. To generate a bubble, the valve was set to divert the liquid flow 
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through the bypass segment, while a wire heater evaporated liquid in the bubble generation 
segment. The bubble volume was varied by adjusting the power to the wire and the heating 
time. Once the desired bubble volume was generated, it was released by rotating the valve 
such that liquid was redirected through the bubble generation segment, pushing the bubble 


into the test section. 


Three-way 
Valve 
Bubble Bypass 
Generation Segment 
Segment 


Heater Wiring ———~ 


Liquid Inlet 


Figure 2.3: Schematic of bubble generation section: flow pattern (left); isometric view with valve 
handle extension (right). 


The bubble rose vertically into a 6 mm ID silicon test section (see Figure 2.4) positioned 
200 mm downstream of the three-way valve where heat transfer measurements and infrared 
flow visualization were made. Pressure taps were located at the inlet and outlet of silicon 
tube so differential and absolute pressures could be measured. The absolute pressure 
transducer (Omega PX209-030A5V) was used to determine the saturation temperature of 


the fluid entering the test section. Immediately after leaving the heated silicon tube, the 
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rising bubble passed through a glass adiabatic section (380 mm downstream of the three- 
way valve) where high speed video was obtained using CMOS video cameras (Phantom 
Miro eX4 and Sentech STC-MBCM200U3V) at frame rates between 500 and 1200 frames 
per second. The high-speed visualization was used to determine the bubble length, study 
the dynamics of the tail, and analyze the bubble shape. The liquid film thickness was 
measured using two techniques: a laser displacement sensor (Keyence LK-G5000) and 
high-speed image analysis. The Keyence sensor was calibrated by measuring a known 
thicknesses of HFE 7100 within a maximum uncertainty of 17 m. Due to movement of 
the sensor relative to the tube when the test apparatus was flown on the aircraft, data from 
the Keyence sensor could not be obtained. The film thickness under these conditions were 
therefore calculated from analysis of the high-speed video captured in the adiabatic section. 
Error in the measurement created by optical distortion was removed by a calibration 
procedure similar to Liu et al. [27] in which a grid of known spacing was placed in the test 
section tube. The apparent spacing was compared to the known spacing to create a 
calibration curve. The uncertainty associated with this technique was estimated to be 25 
pm. Further information on the calibration and uncertainty analysis of the film thickness 


measurement is provided in Appendix A. 
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Figure 2.4: Schematic of heat transfer test section. 


Bubbles were condensed and the liquid subcooled in a counterflow heat exchanger 
where the secondary fluid was cold water. A bellows-type accumulator was included after 
the condenser with the dry side open to the ambient to maintain the system pressure at 
nominally 1 bar in the lab, or the aircraft cabin pressure (typically 0.76-0.83 bar) during 


parabolic flights. Before re-entering the gear pump, the fluid was sent through a de-gassing 
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membrane (Liqui-Cel SuperPhobic). A vacuum pump connected to the membrane for 
several hours prior to data being collected was used to degas the liquid. 

Transducer data was collected using a 24-channel data acquisition system (Omega 
OMB-DAQ-3000) and recorded at a rate of 100 Hz through a LabVIEW™ interface. T- 
type thermocouples were installed at various locations in the flow loop for data analysis, 
calibration, and safety purposes. Uncertainty for the thermocouples was calculated to be + 
0.12°C. Heat transfer measurements and flow visualization were made using an IR camera 
(Electrophysics Silver 660M) at a frame rate of 246 Hz. Typical uncertainties for the major 


instrumentation and some of the reduced parameters are summarized in Table 2.2. 


Table 2.2: Typical uncertainties for important system parameters and measurements. 


Parameter Uncertainty 
G [kg/m?s] 5.2 

Paps [millibar] 1.3 

6 [pm] 17 or 25 
Team [°C] 0.14 

kp [W/m-K] 0.01 

asi [1/m] 6.5 

Qp [1/m] 192 

Tsat [°C] 0.14 


2.2 IR TECHNIQUE 


Heat transfer measurements and flow visualization were obtained using an IR thermometry 
technique developed by Kim et al. [47] that takes advantage of the transparency of silicon 
in the mid-IR range (3-5 pm). HFE 7100 passed through the 6 mm ID (8 mm OD) single 


crystal silicon tube which was doped to allow for resistive heating. The input power to the 
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tube was varied by adjusting the voltage of a high-voltage power supply. The inner wall 
was of the tube was coated with a 57 ym layer of polyimide tape (k, = 0.12 W/m-K) as 
shown in Figure 2.5 One half of the inner circumference was then covered with an IR 
opaque paint containing carbon black (Nazdar GV111), which allowed an effective inner 
wall temperature to be measured through the silicon and polyimide layer. Two strips of the 
painted polyimide tape were also attached to the outer wall of the tube so the outer wall 
temperature could be measured. To complement the heat transfer measurements, the flow 
was visualized using a set of six gold-plated mirrors (Figure 2.5, left image) arranged such 
that flow visualization and heat transfer measurements could be captured using a single 
camera. A representative IR image illustrating the two “halves” of the tube is shown in 
Figure 2.6. The top half was used to visualize the flow within the channel, while the bottom 
half was used for temperature measurements. Note the two white (warmer) strips on the 


bottom half which are the two black tape strips on the outside of the tube. 


8mm OD | 


Flow Visualization 
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Mies Trance Silicon Tube 


Black Coating 


Gold Mirrors 


Figure 2.5: Mirrors to provide simultaneous heat transfer measurements and flow visualization (left 
image), and cross-sectional view of silicon tube with coated polyimide tape (right image). 
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Figure 2.6: Representative IR image showing tube "halves" observed using mirror arrangement. 
The top half was used for visual observation of the internal flow, while the bottom half allowed for 
temperature measurements to be made. 


The temperature profiles within the multilayer were calculated using a 1-D heat 
conduction equation where the governing equations for the silicon (eqn. (6)), adhesive 


(egn. (7)), and polyimide layers (eqn. (8)) are given by 


aT . 
Psilp,si 5, = ksiV?T + Gsi (6) 
aT , 
Padp,ad Ot = KaqV°T (7) 
aT 
Popp op kpV°T (8) 


with boundary conditions T=T;,.(z,t) at the outer silicon tube wall and T=T;,((z,t) at the 
liquid/polyimide interface. A schematic of the multilayer is shown in Figure 2.7 with the 
two boundary conditions labeled along with representations of the thermal gradients 


present within the layers. 


28 
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Figure 2.7: Schematic of silicon tube and tape multilayer showing the thermal gradients present 
due to heating of the tube and cooling from forced convection and boiling. 


Due to absorption, reflection, and emission of energy by the various layers, the observed 
temperature of the inner and outer black surfaces were not indicative of the actual 
temperatures. Instead, a set of optical calculations was coupled with the conduction 
equations to determine the correct boundary conditions. The temperature of the tube outer 
surface was obtained by an iterative process involving two unknowns. The theoretical 
energy collected by the camera while observing the outer tube surface could be calculated 
by summing the energy from each of the sources, 

Eco = 1 — &s) Ew + EsEs.0 (9) 
where €, is the emissivity of the black paint, E.. = Fy,-2, oT, is the energy contributed by 
the surroundings, and E,, = F;,-2, oT;', is the blackbody radiation emitted by the outer 


black surface. F,-,, is the fraction of total blackbody energy collected by the camera over 
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its bandwidth (A, — Az = 3 —5 pm) ata particular temperature. As a result, Fy, 3, varied 
for each source term in eqn. (9) and was determined by integration of Planck’s function 
over the wavelengths and temperature of interest. T.. was determined experimentally by an 
in-situ calibration in which liquid HFE 7100 at a known temperature flowed through the 
tube at a high mass flow rate. Due to the high thermal conductivity of silicon, an effective 
Tz. was calculated such that T;, was equal to the fluid temperature. This procedure was 
repeated for several fluid temperatures to ensure an accurate calibration. 

With the surroundings source term defined, the remaining unknowns were T;, and 
F,,-a,,» which is a function of T;,. A temperature for the tube outer surface was assumed 
and the theoretical energy collected by the camera calculated. After comparing to the actual 
collected energy, T; was updated and the process repeated until the change in temperature 
was acceptable (AT, < 10~° °C). 

Determination of the polyimide/liquid interface temperature was more complex as the 
absorption, emission, and reflection of the multilayer components were accounted for. The 
theoretical energy collected by the camera while observing the inner black surface was 
calculated using the equation, 

Foi = Po-cHo + €si-cEsi + Ep-cE'T + Ts—cEs,i (10) 


where Ex; = i "as; Fy _-a,O1Tsi (x) ]*exp(—as;x)dx is the energy emitted by the silicon 


that reaches the air/silicon interface, Er = i ar Fy,-,,0[Tr(x)]*exp(—arx)dx is the 
energy emitted by the combined adhesive and polyimide layer that reaches the 
adhesive/silicon interface, and E,; = Fy,-a, oT}; is the blackbody radiation emitted by the 
inner black surface. The properties P.—¢, €si_c, Er—c, and T,_, are the effective reflectivity 


of the multilayer, effective emissivity of the silicon, effective emissivity of the combined 
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adhesive/polyimide layer, and effective transmissivity of the multilayer, respectively. Their 
derivations can be found in Kim et al. [47]. A separate value of T.. was determined for the 
inner wall calculations by a calibration technique in which the inner and outer wall 
temperatures were measured while single-phase liquid was passed through the silicon tube 
with constant power applied to the tube. Heat losses from the tube to the surroundings and 
test section housing were approximated using the tube outer wall profile. Further 
information on the loss calculations are provided in the appendix (section A.2.7). The heat 
losses were subtracted from the input power (measured using the four-wire method) to 
determine the net heat absorbed by the fluid. An energy balance using thermocouple 
measurements upstream and downstream of the test section were in agreement with these 
calculations. An effective T.. was then calculated such that the average heat flux along the 
tube length, determined from the measured wall temperatures and MATLAB analysis, 
agreed with the net heat flux value. This procedure calibrated out variations in optical 
properties with temperature, reflections from within the tube, and thermal contributions 
from the surroundings. 

With E,, defined, the theoretical energy calculated by eqn. (10) was compared to the 
actual energy collected by the camera. If a difference was observed, an updated T;,; was 
determined. The inner and outer tube temperatures were used as a boundary condition for 
the conduction equations to determine a new temperature profile. This procedure was 
repeated at each time step until T;; no longer changed by a significant value (AT;, < 
10~° °C). The final calculated temperature profile at each t was used as the initial 


assumption for the following time step. 


31 


The instantaneous heat flux at each pixel was determined by calculating the temperature 


. ea eee - aT 
gradient at the polyimide/liquid interface, or q" = —k, aE . For general two-phase 
X=X5 5 
heat transfer measurements, the heat transfer coefficient was defined as h = Ean, 
si” 4 sat 


where q” and T;,; are defined at each axial pixel location by the calculations above. The 
saturation temperature (7,,;) was found using the test section absolute pressure, according 


to saturation data provided by 3M (manufacturer of HFE 7100). The heat transfer 


coefficient for Taylor bubble experiments was defined as h = where 7; is the 


q 
(Tsi-Tm)’ 
mean fluid temperature at each axial location calculated by an energy balance on the fluid 
using the measured heat flux. This definition was chosen because slightly subcooled 


conditions (1—3°C) were required to prevent undesired nucleation in the pre-heater. 


2.2.1 Validation 

The IR thermography technique was validated through single-phase and two-phase testing 
in a vertical upward flow configuration. For turbulent single-phase flow, liquid passed 
through the heated silicon tube and the heat transfer coefficient along its length was 
measured. The experimentally measured heat transfer coefficient was compared to the 
Dittus-Boelter equation corrected for the thermally developing flow at the beginning of the 
heated silicon tube using the factor proposed by Al-Arabi [48]. The experimental data 
shown in Figure 2.8 is in good agreement with the correlation. Additional validation of the 
technique for the case of laminar single-phase flow was completed by measuring the 
temperature of the tube inner wall at two heat fluxes. Heat losses attributed to natural 
convection on the tube outer wall and axial conduction along the tube into the test section 


end caps were calculated to be 32+4% at these conditions (23+3% to axial conduction, 
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9+1% to natural convection). Losses due to axial conduction were approximated using the 
temperature gradient at the downstream end of the tube, measured using the infrared 
technique. Natural convection losses were estimated using a correlation by Churchill and 
Chu [49] for free convection from a vertical wall, where the surface temperature was taken 
as the average outer wall temperature of the tube. An energy balance on the test section 
conducted by measuring the inlet and outlet liquid temperatures during heating confirmed 


- 3 
oP) to the square of the 


the calculations. The ratio of Grashof number (Gr = 
liquid Reynolds number was of the range 0.68 < = < 1.17, indicating that natural 
l 


convection and forced convection likely contribute to heat transfer from the inner tube wall, 
according to Incropera et al. [50]. In the absence of available correlations from the 
literature, Incropera et al. suggest that the mixed convection Nusselt number (Nuy) be 
computed using the form Nuy = Nug + Nuy, where Nu, is the forced convection 
Nusselt number and Nuy is the natural convection Nusselt number, both computed using 
existing correlations for the specific geometry. In this case, Nuz was approximated from 
Shah and London [51] and Nu, from Davis and Perona [52]. It was noted that the best 
correlation of data is often obtained for n = 3, which was chosen for this case. Accounting 
for heat losses, the predicted and measured inner wall temperature (Ty) offset by the tube 
inlet temperature (Tin) is plotted as a function of axial position in Figure 2.9 showing 
relatively good agreement within the uncertainty of the measurements and predictions. The 
decline in wall temperature at the tube exit is dominated by wall axial conduction losses. 
Due to limitations in microgravity and hypergravity duration, extensive validation of 
the experiment was not able to be conducted under those conditions. However, energy 
balance tests, as described above, were repeated in the laboratory based on the conditions 
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experienced on the aircraft. Notably, the effect of natural convection variation on the outer 
tube wall could not be studied, but two comments in this regard may be made. First, 
according to the correlation of Churchill and Chu [49], the average heat transfer coefficient 
for free convection over a vertical plate is related to gravity by the 1/4" power. As a result, 
an increase in gravity level by a factor of 1.8 (hypergravity) yields an increase in the heat 
transfer coefficient by 16%, yet only an overall increase in power lost by 1.4% (9% x 

1.16 = 10.4%). Second, natural convection loss for a/g = 0.01 was assumed to be 
similar to a/g = 1, given that the data collected under microgravity was taken 2-3 seconds 
after the transition from a/g = 1.8 to a/g = 0.01 and currents associated with free 
convection may have still existed. All data was processed using the corresponding heat loss 
calculations for that gravity level. 

Two-phase validations were conducted by comparing data obtained in the chum and 
annular flow regime to data collected by a group at IMFT in Toulouse, France using a test 
apparatus of different design but operated under similar conditions [53]. Figure 2.10 shows 
that data collected from both experiments are in agreement with each other as well as the 


correlations by Chen [54] and Cioncolini and Thome [55] at low vapor qualities. 
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Figure 2.8: Experimental heat transfer coefficient compared to the Dittus-Boelter correlation with 
the Al-Arabi correction for the parameters: HFE 7100, Re=5545, Ts3=20°C, q’=20 kW/m?. 
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Figure 2.9: Experimental heated wall temperature offset by the liquid inlet temperature as a function 
of axial length along the tube compared to correlations for mixed convection. The conditions were: 
G=50 kg/m?s, Re=790, Tin=45°C, and q"= of 1.4 kW/m? and 2.5 kW/m?. 
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Figure 2.10: Comparison of two-phase data collected by UMD and IMFT with similar experiments 
at conditions - UMD: G=100 kg/m’s, q"=10 kW/m’; IMFT: G=100 kg/m’s, q’=9.8-36 kW/m? 
(figure adapted from Narcy et al. [53]). 


2.2.2 Uncertainty analysis 

The uncertainty associated with heat transfer measurements can be divided into two parts: 
bias error and random error. Bias uncertainty is attributed to possible error associated with 
the determination of material optical properties, thicknesses, thermocouple readings, and 
camera temperature readings. An additional uncertainty relating to the approximation of 
heat losses was also present, but considered as a separate source of bias error. Major 
sources of bias uncertainty and their values are summarized in Table 2.2, including the 
absorptivity of the silicon tube and polyimide tape, thermal conductivity of the polyimide 
tape, and IR camera temperature measurement. The bias uncertainty in measured heat flux 
was calculated to be 6q"=2.3 kW/m? over the range of applied tube heat fluxes in this work. 
Bias uncertainty in the heat transfer coefficient varied linearly with applied heat flux, 
ranging from 6h=0.30 kW/m°K at q”"=2.0 kW/m? to 6h=0.73 kW/m°K at q"=11.7 kW/m’. 
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The values of the bias uncertainty were the same order of magnitude of the measurements 
themselves in most cases, but their effect can be partially negated by an offsetting 
procedure due to the fact that the error sources are fixed values or measurements from set 
calibration curves (with the exception of heat loss calculations). Therefore, if the measured 
heat flux for single-phase flow observed just ahead of the bubble is subtracted from the 
heat flux measured when the bubble is present, then the heat flux enhancement due to the 
presence of the bubble can be defined as q"enn = (q" + (6q")pias) — (Q"sp + 
(Sq"sp)pias) £ (69") neat toss = 9" — A" sp + (6q" — 6q"sp)pias + (69")heat toss: If the 
bias uncertainty is same for both measurements, the heat flux enhancement reduces to 
denn = V' — G"sp £ (6q") neat loss: If the bias uncertainties are unequal, the difference 
must be incorporated into the overall uncertainty by propagating the errors. 

The random uncertainty was calculated by means of propagation of error for averaged 


measurements proposed by Taylor [56] in which the averaged uncertainty decreases with 


respect to the individual measurement uncertainty by a factor of 1/VN, where N is the 
number of averaged measurements. Each temporally averaged bubble heat transfer profile, 


a discussion of which will follow later, is considered an individual measurement whose 


ort 
d 


random uncertainty is characterized by 6, = Ti 


where o,. is the standard deviation of the 


heat flux or heat transfer coefficient about the mean value at each location along the bubble 
and t is the parameter associated with measurements which are assumed to follow the 
Student’s ¢ distribution. Hereafter, error bars presented in heat transfer figures denote the 
random uncertainty associated with the particular measurement as well as the uncertainty 


remaining from bias offsetting of measurements. 
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2.3 DATA COLLECTION 


The procedure for normal gravity data collection began with the initialization of the 
experimental rig and components therein. Once circulation was established within the flow 
loop, heating of the liquid to near saturation was begun and the system was allowed to 
come to equilibrium over a period of roughly 30 minutes. Just prior to testing, a “constant 
temperature” data point was taken in which liquid of known input temperature was passed 
through the test sections at a high flow rate as a means of calibration for the infrared video 
and thermocouples. The silicon test section was then set to the desired applied wall heat 
flux and several bubbles were generated to ensure proper operation of the system. The 
procedure for each data point was as follows: 1) a bubble of desired volume was created in 
the generation section; 2) recording of the test conditions was started via LabVIEW™; 3) 
infrared and film thickness data collections were begun; 4) the bubble was released and 
tracked through the various sections to ensure no secondary bubbles were present; 5) 
recording of test conditions, infrared, and film thickness data was stopped. This process 
was completed for 10-15 bubbles during each collection session. Data was analyzed using 
a series of MATLAB [57] scripts yielding calibrated test conditions, heat transfer 
measurements, and bubble film profiles. 

The experiments completed under a normal gravity environment provide insightful 
information into the heat transfer and flow characteristics of Taylor bubbles rising under 
the conditions available in the laboratory. The bubble drift velocity, Ua, however, is a 
parameter that will be shown to play a role in determining the bubble shape and the heat 
transfer in the following sections, but is not easily varied. A slight manipulation of egn. (2) 


yields a relation for Ug, which is seen to be dependent on the parameters C, U}, and Up, 
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Ug = U, — U, = (C — 1)U, + Up 0 (11) 
where Up, has been found to be function of the difference in the liquid and vapor densities, 
gravity level, liquid viscosity, liquid surface tension, and tube diameter. The range of 
bubble drift velocities available for the laboratory experiments was limited as the dominant 
term in egn. (11) is Ub. Modification of the test section tube diameter or exchanging the 
test fluid would provide an expansion of the drift velocity range, but would prevent direct 
comparison of old and new data sets. Varying the gravitational acceleration experienced 
by the bubble delivers a less obvious, yet more straightforward solution. 

Several means exist to increase or decrease the gravity level in the experimental 
environment. Drop tower tests are relatively inexpensive and can be repeated several times 
per day, but the microgravity duration is too short (2-3 seconds) to allow for a bubble to be 
released and tracked through the test sections. Conducting experiments on a sounding 
rocket or the International Space Station, while providing longer microgravity periods 
(minutes to months), are cost and size prohibitive. A median solution is a parabolic flight 
campaign in which researchers are provided 20-30 seconds of microgravity and 
hypergravity, while being able to work directly with their test rig. Experiments were, 
therefore, conducted aboard NASA’s C-9 aircraft (shown in Figure 2.11a) which was based 
at Johnson Space Center in Houston, TX. The campaign consisted of four flights with each 
consisting of forty parabolic trajectories similar to the one shown in Figure 2.11b. Each 
parabola provided up to 22 seconds of microgravity (0.02+0.01g) and up to 30 seconds of 


hypergravity (1.8+0.01g) during which data could be collected. 
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Figure 2.11: Description of parabolic flight campaign: a) NASA C9 aircraft; b) parabolic flight 
profile (photo and graph courtesy of NASA). 


The data collection procedure during parabolic flights was slightly abbreviated 
compared to experiments conducted in the laboratory. Due to aircraft fuel and flight plan 
restrictions, the available time to initialize the test rig and preheat the flow loop was 
reduced to 10-12 minutes. The reduced pressure within the cabin, however, provided a 
lower required heating temperature as the saturation temperature of the fluid was 
decreased. As a result, the restricted warming time was still sufficient to allow the system 
to reach an equilibrium temperature. Individual data points were collected using the same 
series of steps described for normal gravity testing for a total of approximately 40 bubbles 
during each flight. The extended duration of testing required “constant temperature” data 
to be collected periodically to account for heating of the infrared camera housing and the 
rig containment volume. Data was later analyzed as described previously. 

Normal gravity and variable gravity results were combined into a complete data set 


which spans a wide parametric space and provides novel insight into heat transfer 
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mechanisms of rising elongated bubbles. These results are presented and discussed in the 


following chapter. 
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CHAPTER 3: RESULTS AND DISCUSSION 


Data was collected in the laboratory and during a parabolic flight campaign which 
consisted of 164 parabolic maneuvers. Each parabola provided approximately 20 seconds 
of microgravity (+0.01g) and 30 seconds of hypergravity (1.8g+0.1g) where g is 
gravitational acceleration. Background liquid velocities within the test section were of the 


range 35 mm/s < U;< 140 mm/s (50 kg/m?s < G < 200 kg/m’s), which correspond to liquid 
Reynolds numbers (Re, = -) between 790 and 3090. The variation of these parameters 
l 


allowed for heat transfer measurements and flow visualization to be obtained at 0.49 < Bo 
< 87, thereby spanning both capillary and Taylor bubble regimes. A summary of the test 


conditions is given in Table 3.1. 


Table 3.1: Summary of data parameters. 


Ua (mm/s) _ G (kg/m’s) Re; a/g Bo Urs(mm/s) Up (mm/s) 

20 100 1580 0.01 0.48 — 90 

41 200 3090 0.01 0.48 - 180 
94 200 3090 0.34 16.5 - 235 
105 50 790 1 48.4 851 142 
124 100 1580 1 48.4 1030 196 
144 50 790 1.8 87.2 1151 180 
163 200 3090 1 48.4 1402 305 
173 100 1580 1.8 87.2 1389 241 
208 200 3090 1.8 87.2 1754 349 
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3.1 BUBBLE DYNAMICS 


Elongated bubbles released from the bubble generation section rose into the silicon test 
section where the heat transfer was measured. High-speed video was then obtained as the 
bubble passed through the adiabatic glass section. The bubble velocity was determined by 
tracking the bubble nose as a function of time. A series of 5—20 bubbles were analyzed and 
averaged to obtain the bubble velocity at each liquid velocity for four gravity levels as 
shown in Figure 3.1. For a given liquid velocity, the bubble velocity was found to increase 
with increasing gravity due to the larger buoyancy force acting on the bubble. A line was 
fit to the velocity measurements at three gravity levels (a/g=0.01, 1, and 1.8) to compare 
the data with correlations of the form proposed by Nicklin et al. [16]. Results for Up,o are 
compared to the predictions of White and Beardmore [18], Rattner and Garimella [20], 
Brown [19], and Viana et al. [13] in Table 3.2. Very good agreement is seen with Brown, 
Viana et al., and White and Beardmore for normal gravity conditions (1g), while the 
correlations tend to slightly under-predict Up,o at hypergravity. As expected, Up,o was very 
small under microgravity (+0.01g) conditions due to the small residual buoyancy force. 


Similar observations were made by Colin et al. [58] for Taylor flow in microgravity. 
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Figure 3.1: Bubble velocity as a function of liquid velocity for a/g=0.01, 0.34, 1, and 1.8. 


Table 3.2: Comparison of measured Up, to those calculated from published correlations. 


Present White and Viana et al. Rattner and 
a/g Study Beardmore [18] Brown [19] [13] Garimella [20] 
0.01 2.9+5.8 0.0 — — 0.0 
1 81.1+4.1 79.6 82.0 82.4 72.4 
1.8 130.2+8.8 112.3 110.4 110.6 108.8 


Calculated values for C at 1g and 1.8g (Bo=49 and 87) were found to fall within the 
classical bounds of 2 and 1.2 for the laminar and turbulent liquid flow regimes, 
respectively. Under microgravity conditions (Bo=0.49), C was measured to be smaller than 
at normal and hypergravity. The variation is likely due to a transition from Taylor flow to 
capillary flow given that the bubble Bo in microgravity lies below the critical values 
suggested by previous authors [18,21,36—38]. Visual evidence of a transition is provided 
in Figure 3.2 where high-speed images of bubbles rising at various gravity levels and 


bubble drift velocities (Ug = U, — U,) are shown. Under microgravity conditions (Figure 
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3.2a,b), the bubbles exhibited rounded tail profiles whose radius of curvature increased 
with drift velocity. Stationary, small amplitude interfacial waves were observed at Ug=20 
mm/s (Figure 3.2a). When the drift velocity increased to 41 mm/s (Figure 3.2b), the waves 
became more frequent, larger in amplitude, and oscillated in the streamwise direction. 
Above Ug=94 mm/s (Figure 3.2c—e), the Taylor bubbles possessed flat tails that wobbled 
as the bubbles rose. The magnitude of the wobble was found to increase with bubble length. 
Figure 3.3 illustrates the tail shape of two bubbles rising in a co-current flow at Ui=37 mm/s 
and a/g=1. The shorter bubble (L,/D=2.1) experienced mild oscillations and the tail profile 
remained relatively perpendicular to the flow direction, while the longer bubble (Lp/D=8.8) 
oscillated significantly and exhibited irregular tail shapes. This trend is in agreement with 
Polonsky et al. [22], who showed that a small increase in oscillation frequency and large 
increase in amplitude accompanied bubbles with longer lengths. 

A more detailed bubble profile analysis was conducted using the open source image 
processing program ImageJ [59]. The nose profiles of bubbles with drift velocities of 
Ug=20, 41, 94, 175, and 214 mm/s are compared in Figure 3.4. The bubbles exhibit 
narrower, more rounded nose contours with increasing Ug. The theoretical nose shape 
suggested by Dumitrescu [60] for a Taylor bubble rising in a stagnant column of fluid 
without viscosity or surface tension at normal gravity is also included. With a 
corresponding drift velocity of approximately 81 mm/s (the value of Uso at 1g), 
Dumitrescu’s profile is seen to fall between the experimental curves for Ug = 41 and 94 
mm/s near the top of the bubble. As the liquid film develops, the prediction moves within 


the measured profile for Ug=94 mm/s but remains within the uncertainty of the experiment. 
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a/g 


a) b) c) d) e) 


Figure 3.2: Images of representative bubbles at varying accelerations and drift velocities: a) a/g = 
0.02, Ug = 20 mm/s, Bo = 0.97, L=17 mm; b) a/g = 0.02, Ua = 41 mm/s, Bo = 0.97, L=18 mm; c) 
a/g = 0.34, Ua = 94 mm/s, Bo = 16.6, L=13 mm; d) a/g = 1, Ua = 106 mm/s, Bo = 48.7, L=13 mm; 
e) a/g = 1.8, Ua = 214 mm/s, Bo = 87.7, L=13 mm. 


46 


z/D 


0 1 
-0.5 0 : 
r/D r/D 


Figure 3.3: Representative tail profiles for two bubbles rising in co-current flow with U=37 mm/s. 
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Figure 3.4: Bubble nose profiles at Ua=20, 41, 94, 175, and 214 mm/s. The theoretical profile of 
Dumitrescu for Ua=81 mm/s is also shown. 
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As was mentioned, small waves were observed for capillary bubbles (Ug=20 and 41 
mm/s) at the vapor/liquid interface near the bubble tail. Profiles for two representative 
bubbles at these drift velocities are shown in Figure 3.5. The wavelength (A) of the waves 
was found to be slightly dependent upon Ug with values of 2.13 + 0.04 mm and 1.79 + 0.04 
mm for Ug=20 and 41 mm/s, respectively. Ratulowski and Chang [61] numerically 


predicted the shape of elongated bubbles as a function of the Capillary number 
(é a= ae), A series of relations were developed for the characteristic wavelength of 


capillary waves at the bubble interface. For Ca=0.0034 (Ug=20 mm/s) and Ca=0.0069 
(Ua=41 mm/s), Ratulowski and Chang predict A to be 2.15 mm and 2.83 mm, quite similar 
to the current experimental results. Edvinsson and Irandoust [62] compared their numerical 
results for 2 to those of Ratulowski and Chang, and found that 1/D was overpredicted by 
the original results for Ca>0.005. Edvinsson and Irandoust’s results were applied to the 
current study yielding a wavelength of 2.41 mm for bubbles with Ua=41 mm/s, which is 
somewhat closer to the current measurement of 1.79 + 0.04 mm than predictions from 


Ratulowski and Chang. 
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Figure 3.5: Bubble shape profiles for Uz=20 mm/s (left) and Ug=41 mm/s (right) under microgravity 
conditions. 


The liquid film thickness for Taylor bubbles (Ug > 94 mm/s) was obtained using either 
the Keyence sensor or high-speed image analysis as described previously. Figure 3.6 
compares data collected in 1g (Keyence sensor) and 1.8g (image analysis) with their 
respective uncertainties labeled on several points. The tail film thickness for bubbles with 
Ug=105, 124, 163, 173, and 208 mm/s are plotted as a function of the bubble length. As 
expected, 6 decreases with bubble length until approximately Lp/D=6 where the thickness 
reaches the fully developed value. This development length (Z) was compared to 
predictions of Campos and Carvalho [23], which range from 13.8 to 16.1 for the current 
experimental results depending on the bubble velocity. Nogueira et al. [24] found that at 


film Reynolds numbers (Rey oo (U, — U,)D/4v) above approximately 80, the Campos 


and Carvalho correlation over-predicted their experimental results by as much as 30%. The 
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current data shown in Figure 3.6 suggests Lp/D~6 is required for the film thickness to 
stabilize, which is not inconsistent with either of these predictions. 

The data is compared to predictions by Brown [19] for fully developed film thickness 
(6s) by two means. First, dg was calculated for each value of Ua using the measured Us, Ui, 
and an assumed acceleration of 1g. This range of 6z is illustrated by a blue area in Figure 
3.6. Also shown is the prediction for Uaz=163 mm/s and a/g=1, which agrees well with the 
measurements at these conditions. The calculations were repeated for Ug=173 and 208 
mm/s at hypergravity, with the results shown by the red area and dash-dot line (Ug=173 
mm/s). It is clear that the relation by Brown over-predicts the film thickness for higher 
gravity levels. In fact, the measurements at a/g=1.8 vary little from those collected at a/g=1. 
Using the correlation of Llewellin et al. [26] for fully developed film thickness of Taylor 
bubbles in stagnant flow, 6 was calculated to be 271 ym for a/g=1 and 264 pm for a/g=1.8. 
Both of these values fall within the 1g prediction range of Brown and illustrate relatively 
little change in film thickness with increasing gravity level, which is consistent with the 


measurements. 
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Figure 3.6: Film thickness at the bubble tail as a function of bubble length at various drift 
velocities for Taylor bubbles (Ug>94 mm/s). 


3.2 HEAT TRANSFER 


The wall heat transfer was measured for each bubble as it passed through the silicon test 
section where steady-state, thermally developing flow was established to ensure 
reproducible boundary conditions. Measurements were offset by the single-phase heat 
transfer observed just before the bubble enters the test section. This means of 
characterization was chosen to emphasize the effect of heat transfer enhancement created 
by the bubble passage, rather than the effect of the background flow conditions. Two heat 
transfer coefficient profiles (h — hsp) obtained at a Reynolds number of Re;=1580 
(a/g=0.01 and 1) are plotted in Figure 3.7 as a function of axial location along the bubble, 
where z/D=0 indicates the bubble tail. The profiles shown are an average of the 
instantaneous profiles measured in the frame of reference of a viewer moving with the 
bubble tail. For both gravity levels, the heat transfer coefficient ahead of the bubble 


ol 


corresponds to single-phase flow. As the liquid thins during bubble passage, a slight 
enhancement to the heat transfer occurs due to increased conduction and convection. 
Waves in the liquid-vapor interface create oscillations in the capillary bubble heat transfer 
(Figure 3.7a) as the fluid accelerates and decelerates within the film. At the capillary bubble 
tail (~z/D=0), a significant decline in the heat transfer enhancement was observed and may 
be attributed to a local thickening of the thermal boundary layer caused by a negative radial 
velocity in this region. Magnini et al. [41] also observed these flow characteristics in their 
numerical simulations of capillary flow within a microchannel. No such oscillations are 
seen in the Taylor bubble profile (Figure 3.7b) where the interface remains smooth. At the 
tail of the Taylor bubble, vortices are generated when the downward moving film plunges 
into the trailing liquid slug, inducing turbulence and a large spike in heat transfer. Vortices 
were not observed in the wake of the capillary bubble and the heat transfer coefficient was 
seen to be similar to single-phase flow. Representative IR flow visualization images of the 
two bubbles are shown in Figure 3.8. The wake (to the left of the image) of the capillary 
bubble was seen to maintain an essentially uniform temperature, suggesting that the near- 
wake region remained laminar and no mixing occurred. Just behind the Taylor bubble, 
however, clear temperature fluctuations (vortices) can be seen, indicating mixing of 


warmer and cooler fluid. 
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Figure 3.7: Representative heat transfer coefficient profiles: a) capillary bubble (Uga=20 mm/s, 
U=70 mm/s, L»=16 mm, a/g=0.01, q”=800 W/m?, hsp=101 W/m7k), b) Taylor bubble (Ug=124 
mm/s, U=70 mm/s, L)=14 mm, a/g=1, q”"=1440 W/m’, hsp=148 W/m’K). 
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Figure 3.8: Representative IR images of the near wake region for bubbles presented in Figure 3.7: 
a) capillary bubble and b) Taylor bubble. 


For Taylor bubbles, the bubble length was found to have little effect on the shape and 
magnitude of the wake heat transfer enhancement as shown in Figure 3.9. The enhancement 
due to five bubbles of different lengths at Uz=163 mm/s differed by an amount similar to 
the uncertainty in the measurements and no visible trend in the profiles was seen. If it is 
assumed that the relative velocity between the liquid film and trailing liquid slug is the 
parameter that dictates wake heat transfer, then this result is not surprising. As was seen in 
Figure 3.6, the liquid film thickness is essentially invariant for Lp/D>4. The bubbles 
presented in Figure 3.9 were very near or above this threshold, suggesting that their film 
velocities and subsequently wake heat transfer should also be similar. Babin et al. [34], on 
the contrary, found that the magnitude and length of the wake enhancement increased with 
bubble length for both laminar and turbulent flows. They attributed the trends to two 
factors. First, the higher film velocities associated with longer bubbles resulted in vortices 
that persisted farther downstream. Second, the increased wobbling of the tail for longer 


bubbles may have resulted in increased wake enhancement. 
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An additional explanation arises if the bubble development length is considered. The 
bubble length required for a fully developed liquid film was calculated for the conditions 
of Babin et al. [34] to be Lp/D=22 using Campos and Carvalho [23] where the film 
thickness was approximated using Brown [19]. The bubble lengths considered by Babin et 
al. ranged from L»/D=2-10, shorter than the Lp/D for fully developed flow, even if Z is 
over-predicted by the correlation. It is possible that the enhancement in wake heat transfer 
observed with increasing bubble length is the result of higher film velocities. As the 
bubbles in the current experiment were essentially fully developed, little variation in the 


wake heat transfer profiles are seen for the results. 
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Figure 3.9: Scaled heat transfer coefficient as a function of distance behind the bubble tail for 
bubble drift velocity Us =163 mm/s, tube heat flux of q"=1.0 kW/m’, average single-phase heat 
transfer coefficient hsp=192 W/m’K, and bubble lengths L,/D=3.2, 4.2, 5, 5.8, and 7.7. 


A 2-D plot of the temporal and spatial wake heat transfer can be used to infer the local 


effect of vortices generated by Taylor bubbles. An example of such a plot and the wake 
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characteristics are shown in Figure 3.10. The resolution of the data was limited by the 
spatial (114 pm, 0.02D) and temporal (4 ms) resolution of the camera, thus creating the 
pattern of colored rectangles. Streaks of elevated heat transfer (labeled by black lines) were 
likely caused by the departure of vortices from the tail. The distance behind the tail at which 
the streaks appear (the penetration length from Kawaji et al. [63], Lp) is seen to vary with 
respect to time, a behavior that can be attributed to two factors as illustrated in Figure 3.11. 
First, the location of the tail (z/D=0) is determined from the average bubble velocity and 
does not account for wobbling of the tail. As a result, the actual location of the tail at each 
time may be slightly ahead or behind the average position, leading to a variation in the 
perceived penetration length. The second factor to be considered is the instability inherent 
in the plunging jet which would likely lead to oscillations in Lp even if the tail position was 


fixed. 
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Figure 3.10: Representative wake heat transfer signature characteristics identification. 
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Figure 3.11: Schematic showing high velocity wall jets (red) plunging into the wake. Separation of 
the jet occurs at the penetration length. Bubble photograph courtesy of Dr. Mirco Magnini. 


Heat transfer coefficient contours for flows with Re=790, 1570, and 3090 at a/g=1 and 
a/g=1.8 are shown in the left images of Figure 3.12 and Figure 3.13, respectively. Plots on 
the right of each figure illustrate the heat transfer coefficient profile at four time steps from 
the corresponding contour on the left. These individual profiles confirm the variation in Lp, 
which is seen to range from 0.5-1.5 diameters for bubbles observed at a/g=1 and 0.7-2.0 
diameters for a/g=1.8. Similar results were obtained by Kawaji et al. [29] for air bubbles 
rising in stagnant kerosene, where Lp ranged from 0.85-1.4 diameters depending on the 
bubble length. Shemer et al. [30,64], while not explicitly describing the penetration length, 
found that the highest time-averaged radial velocity in the wake occurred at approximately 
1 diameter from the tail, also consistent with the current results. 

An increase in average Lp was observed for increasing Re; (and subsequently Uf) at each 
gravity level. This behavior was potentially the result of the increased momentum 


possessed by the film as it plunged into the wake. Results of Babin at al. [34] exhibited 
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similar trends, but an explanation was not discussed. Kawaji et al. [29] suggested 
(referencing Lamb [65]) that increasing the relative velocity between the film and the wake, 
as occurred with increasing Re, would increase the rate of instability growth and 
potentially cause shorter values of Lp. This theory, however, based on Helmholtz instability 
analysis, does not account for the increased velocity of the bubble with Re;, which may 
result in larger values of Lp as the bubble moves away from the injection point, despite the 


increased growth rate of the instability. 
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c) Ug=163 mm/s, g"=1700 W/m?, hsp=192 W/m?K 
Figure 3.12: Wake heat transfer coefficient contours (left) and representative profiles (right) at 
a/g=1 for various conditions. 
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c) Ug=208 mm/s, g”=1320 W/m”, hsp=146 W/m?K 


Figure 3.13: Wake heat transfer coefficient contours (left) and representative profiles (right) at 
a/g=1.8 for various conditions. 
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The vortex frequency, fy, was determined by counting the number of streaks in each 
contour and dividing by the time over which they occurred. At each Up, fy was found by 
averaging the calculated frequencies for a number of contour plots (between 2 and 13 plots) 
at each flow condition. Attempts were made at performing Fourier analyses on the streaks. 
However, due the limited length of the acquired signals (resulting from the limited test 
section length) and quasi-periodicity of vortex generation, the results were inconclusive. 
The calculated frequency results should serve as an approximation of the wake turbulence. 
The frequency increased monotonically with increasing Up, as seen in Figure 3.14, and was 
independent of gravity level for a/g=1 and 1.8. It is possible that the two bubbles observed 
under Martian gravity (a/g=0.34) fall within a transitional region between capillary and 
Taylor bubble hydrodynamics and therefore do not align with the higher gravity data. A 
linear fit to the current a/g=1 and 1.8 data suggests that a critical plunging velocity (Up,cr) 
of 127 mm/s exists where vortices begin to appear. For this experiment, Up,:r=127 mm/s 
corresponds to a critical drift velocity of Udacr=34 mm/s and a critical drift Reynolds 
number of Rey, .. = 741. Pinto et al. [28] proposed Rey, ,. = 525 as the boundary for 
incipience of vortex shedding, which compares fairly well with the current result. The 


frequency can be generalized using the Strouhal number (St), defined here as St = 


f6 


—_——., to be St=0.018+0.003. 
(Up—Up,cr) 
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Figure 3.14: Wake vortex frequency as a function of the relative velocity between the film and bulk 
liquid velocity for a/g=0.34, 1, and 1.8. 


The velocity with which the vortices moved away from the tail in the bubble reference 
frame, U,,, was determined from the streaks in Figure 3.12 and Figure 3.13 by measuring 
their slope in the near-wake region. The length of the near-wake region varied depending 
on the liquid velocity and gravity level. At Ug=105 mm/s (Figure 3.12a), for example, the 
near-wake region was defined as approximately —1.2 < z/D < —0.5, while at Ug=208 
mm/s (Figure 3.13c) it was —2.8 < z/D < —1.5. U;j is plotted in Figure 3.15 as a function 
of the wall velocity Uy,. Vortices are shown to move faster than the wall and increase 
monotonically regardless of gravity level. This is not surprising given that the liquid 
velocity near the wall is dominated by the remnants of the plunging jet, as shown in Figure 
1.4. The fluid motion responsible for the streaks likely remains close to the wall and is 


convected downstream at a velocity greater than the wall velocity. 
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Figure 3.15: Velocity of vortices as they move away from the bubble tail (fixed reference frame at 
bubble tail) with respect to bubble velocity at a/g=1.8 and a/g=1. 


The various time-resolved wake characteristics that have been gleaned from Figure 3.12 
and Figure 3.13 can be used to describe time-averaged heat transfer results. To further 
compare the effect of drift velocity and gravity on the wake heat transfer, data were 
averaged with time for 5-20 bubbles at each test condition to obtain the profiles shown in 
Figure 3.16. The near-wake region (Figure 3.16a) is characterized by the rise of the heat 
transfer coefficient enhancement behind the bubble tail. The penetration length is seen to 
increase with Ug for both gravity levels, with a/g=1.8 exhibiting slightly longer Lp than 
a/g=1. Both observations may be attributed to the increase in liquid film momentum with 
increasing drift velocity and gravity level. 

The peak heat transfer enhancement is also seen to vary with drift velocity, generally 
decreasing with increasing Ug at each acceleration. It can be noticed that the peak behavior 
over both accelerations appears to coincide with the relative distribution of Lp at each drift 


velocity. At a/g=1.8, Lp is seen to vary widely (Figure 3.13), thereby creating more diffuse 
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and smaller magnitude peaks than a/g=1 (Figure 3.12). Definitive elucidation of these 
trends requires coupling of local flow field measurements, which were not available in this 
study. 

The far-wake region (Figure 3.16b) shows a continuation in enhancement decay to the 
single-phase value far downstream of the bubble. The distance behind the tail at which 
single-phase heat transfer is re-established (the calming length, L-) is seen to depend on the 
liquid Reynolds number and the gravity level. This behavior may be attributed to the 
redevelopment of the thermal boundary layer in the wake after its disruption by vortices as 
was shown by Babin et al. [34]. For both 1g and 1.8g accelerations, L-/D was smallest at 
Re;~800, increased to its largest values at Re=1570, then decreased to moderate lengths at 
Re}=3058. Assuming laminar flow, increasing Re; from approximately 800 to 1570 would 
have resulted in longer thermal and hydraulic redevelopment lengths. Increasing Re; further 
to 3080, where the flow was likely transitional, would have caused a decrease in the 
redevelopment length due to the onset of turbulence. The growth in L/D with acceleration 
may be attributed to the increased influence of natural convection at hypergravity 
(Rapisg/Rapig © 2.8, 2.1, and 2.3 for Re, = 800, 1570, and 3090, respectively) which 
effectively increased the Reynolds number near the wall and led to a longer development 


length. 
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Figure 3.16: Heat transfer coefficient enhancement as a function of distance behind the bubble tail 
for a) near-wake region and b) far-wake regions at the conditions: Ug=105, 124, 163, 144, 173, and 
208 mm/s; tube heat flux of q"=1400, 1470, 1700, 1140, 1330, and 1320 W/m’, respectively; 
average single-phase heat transfer coefficient of hsp=130, 148, 192, 104, 112, and 146 W/m°’K, 
respectively. 
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3.3 COMPARISON TO NUMERICAL SIMULATIONS 


The experimental study presented above succeeded in determining major characteristics of 
the heat transfer and flow patterns around rising elongated bubbles. There were, however, 
limitations to the experiment that prevented portions of the fundamental physics from being 
fully understood. Specifically, the local, time-resolved flow and temperature fields in the 
near-wake region, which are invaluable to characterizing the movement of vortices, were 
unobtainable given the experimental setup. Results such as these are more easily captured 
with the help of numerical simulations from which many parameters may be detailed. For 
this reason, a collaboration was formed with Dr. Mirco Magnini working with Prof. John 
Thome at Ecole Polytechnique Fédérale de Lausanne (EPFL) in Lausanne, Switzerland. 
Dr. Magnini has previously developed simulations in ANSYS Fluent for slug flow in 
microchannels that have compared well with the literature [41,42,66]. A detailed 
explanation of the simulation numerical framework can be found in Dr. Magnini’s 
microchannel publications and in a recent UMD/EPFL collaborative conference paper [67]. 
His numerical framework was adapted for use with the current experimental setup to 


provide a comparison to the collected data and to aid in understanding its meaning. 


3.3.1 Simulation Parameters 

In an effort to match numerical and experimental results as closely as possible, a set of 
flow, heating, and bubble parameters were agreed upon to serve as boundary and initial 
conditions for the simulations. These details were determined from the experiments 
completed as part of this study and are summarized below in Table 3.3. The numerical 
boundary conditions were as follows. Laminar fully-developed (Poiseuille) hydraulic 


conditions were assumed at the tube entrance for both diabatic and adiabatic simulations. 
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Bubbles of specified length were initialized at the entrance and simulated to rise through 
the channel length (16D for adiabatic, 40D for diabatic). For diabatic simulations, the tube 
wall was heated with a constant heat flux. Natural convection was enabled to account for 


buoyancy effects in the liquid. 


Table 3.3: Summary of comparative experimental and numerical parameters. 


Parameter Value 
G [kg/m?s] 50 and 100 
q”" [W/m?] 0 and 1380 
Lp [mm] 18 and 36 
g [m/s?] 1 and 1.8 


3.3.2 Shape and Hydrodynamics Comparison 

As was illustrated in the experimental results discussed previously, the behavior of the 
bubble and its surrounding flow field heavily influence the bubble’s effect on the heated 
tube through which it rises. Consequently, it is important that the bubble shape and 
hydrodynamics captured by the simulations closely agree with experiments so that heat 
transfer comparisons can be appropriately made. Representative bubble images from each 
case are compared in Figure 3.17 for G=50 kg/m?s at normal gravity. The experimental 
bubble was visualized in the adiabatic section, while the bubble from the simulation was 
observed in a heated tube. The heat flux applied to the wall was small (g"=1380 W/m?) and 
was assumed to not affect the bubble shape. It is clear that the general shapes of the bubbles 
are quite similar, including the deformation of the tails as they wobbled. More quantitative 


comparisons were also made, including the development of the liquid film. The 
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experimentally measured film thickness is plotted as a function of distance from the bubble 
nose (the measured bubble is shown above the plot) in Figure 3.18 along with the 
corresponding thickness from the simulations. Despite the difference in heating conditions, 
the film thicknesses are seen to agree within the uncertainty of the experimental 


measurements. 


a/g 


Flow 


a) b) 
Figure 3.17: Comparison of representative bubbles under normal gravity conditions with Uj=37 


mm/s from: a) experiments (L=17 mm, adiabatic); b) simulations (L=18 mm, q"=1380 W/m’). 
Simulation image courtesy of Dr. Mirco Magnini. 
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Figure 3.18: Comparison of measured and simulated film thickness for a short bubble (L-.,=14 mm, 
Lnum=17 mm) rising within liquid with U=37 mm/s at a/g=1. The bubble from which the 
experimental film thickness was measured is shown above the figure. Numerical data courtesy of 
Dr. Mirco Magnini. 


The bubble velocity was compared for various liquid mass fluxes at both adiabatic and 
diabatic tube heating conditions. Bubble tracking data from the diabatic simulations 
illustrated a decrease in velocity as the bubble rose in the tube due to flattening of the liquid 
velocity profile caused by natural convection. A variation was not seen for experimentally 
measured velocity, but the distance over which the bubble was viewed was shorter than the 
simulations, which may have truncated the decelerating behavior in the data. Due to this 
discrepancy, a representative value of Up was determined for the simulations by averaging 


the velocity profile over the length that would be visible had it been measured using the 
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experimental apparatus. It can be seen in Figure 3.19 that regardless of the thermal 
boundary condition, both simulation and experimental velocity values agree well with each 


other and linearly increase with U), as expected. 


400 |: @ Experimental (adiabatic), a/g=1 | 
© Experimental (g"=1380 W/m), a/g=1 
@ Numerical (adiabatic), a/g=1 
O Numerical (g”=1380 Wim’), a/g=1 J 
—~ 300} ee 
@ =. 
=, 200} = : 
) 
S ae 
a: 
100.-- 
0 1 it 1 
0 50 100 150 


U, (mm/s) 


Figure 3.19: Comparison of experimentally measured bubble velocities to those obtained from 
numerical simulations. Numerical data courtesy of Dr. Mirco Magnini. 


3.3.3 Heat Transfer Comparison 

The general agreement in bubble hydrodynamics between the two methods of study 
provided encouragement that the simulation results were representative of the experiments. 
Further evidence was attained by comparing wake heat transfer contours in Figure 3.20 for 
G=50 kg/m’s, a/g=1 and q"=1380 W/m’. It can be seen that, like the experiments, the 
contours attained from the simulations were characterized by streak patterns where the 
highest heat transfer coefficient occurred near the tail. The streaks exhibit a higher velocity 


nearer the bubble (near-wake), then transition to a lower velocity downstream (far-wake). 
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Calculated vortex velocities from the simulations were found to agree well with those from 
the experiments (near-wake: Uy), ¢x)=224+11 mm/s, Uj .jm=210+13 mm/s; far-wake: 


Uy exp=119+6 mm/s, Uj, 5¢m=115+16 mm/s). The penetration lengths were also found to be 


similar, ranging from approximately z/D = —0.5 to z/D = —1.5 in both cases. 
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Figure 3.20: Comparison of wake heat transfer coefficient profiles obtained from experiments and 
numerical simulations for G=50 kg/m?s, q”=1380 W/m?, and a/g=1. Numerical data courtesy of 
Dr. Mirco Magnini. 
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Wake heat transfer comparisons made using the contours, while serving to validate both 
experimental and numerical results, do not take full advantage of the simulation data. To 
more completely understand the physics of the bubble wake, the instantaneous flow and 
temperature field can be plotted from the simulations, as shown in Figure 3.21 where 
z/D = 0 indicates the tail position and z/D < 0 is the wake. The velocity field is indicated 
by the arrow vectors and the temperature field is illustrated using the color map. 
Temperatures were offset by the inlet temperature of the heated region, where the liquid 
was subcooled by approximately 4°C. Near the top of the figure (—0.5 < z/D <0), a 
warm, circular jet from the film can be seen to plunge into the cooler liquid behind the 
bubble. Note that as the contour was created from a cross-sectional view, the circular jet 
appears as two jets, one on each side of the tube (z/D = +0.5). The jet becomes unstable 
as it moves downstream and is sheared from the wall, allowing colder fluid to replace it. 
Gradually the liquid behind the bubble mixes to create a more uniform temperature profile 
within the tube. The effect of the temperature and flow evolution was evaluated in Figure 
3.22, where the local heat transfer coefficient at the tube wall is plotted vs. axial distance 
for the two “halves” of the channel. 

The heat transfer coefficient is observed to be relatively small while the warm jet 
remains intact near the tail of the bubble. The jet acts as a barrier for the wall by preventing 
cool liquid from nearing the surface, thereby maintaining a thick thermal boundary layer. 
It is only when the jet begins to break up, at approximately z/D = —0.75 for both walls, 
that mixing of the fluid occurs and the heat transfer rises. Peaks in the heat transfer 
coefficient coincide with regions where the fluid motion transports pockets of cool liquid 


to the wall creating a high temperature gradient. 
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Relating this insight to experimental results, the heat transfer behavior in the wake 
contours can be better explained. First, the penetration length was directly related to the 
wobble of the bubble tail and the destabilization of the circular jet. Due to the complex and 
transient nature of the hydrodynamics, it is not surprising that Lp varied with time in Figure 
3.12 and Figure 3.13. Second, the distribution of heat transfer coefficient with respect to 
z/D can be explained by the evolution of the temperature field within the tube. Just after 
the tail, where the highest heat transfer was observed, the liquid that was brought to the 
wall from the center of the channel was relatively cold compared to the liquid near the tube 
surface. Farther downstream, mixing created a more uniform temperature field and a 
smaller temperature difference existed between the wall and the bulk fluid. This, as well as 


the decay of turbulence, led to a steady decline in the heat transfer coefficient. 
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Figure 3.21: Instantaneous flow and temperature field simulated for a bubble with the conditions: 
G=50 kg/m’s, L/D=6, q’=1380 W/m2, Tin=55.2°C, Tsa=59.5°C, and a/g=1. Data courtesy of Dr. 
Mirco Magnini. 
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Figure 3.22: Local heat transfer coefficient plotted as a function of axial position in the wake for 
the a) left and b) right sides of the cross-section. The averaged bubble parameters were G=50 
kg/m?s, L/D=6, q"=1380 W/m’, Tin=55.2°C, and a/g=1. Data courtesy of Dr. Mirco Magnini. 
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In conclusion, the comparisons presented above showed that experimental results can 
be well simulated by Dr. Magnini’s numerical framework. Agreement was seen in both 
heat transfer and hydrodynamic characteristics for the simulated conditions currently 
available. Furthermore, local wake temperature data from the simulations helped to explain 
the behavior seen in the heat transfer contour plots. Additional simulation cases will be 
completed as part of the collaboration to determine the effects of liquid mass flux and 


gravity level on the heat transfer around a single bubble. 


3.4 IMPLICATIONS FOR SLUG FLOW MODEL DEVELOPMENT 


As one of the main goals of this research was to guide the development of more realistic 
and accurate models for slug flow heat transfer, it is important to discuss how the results 
contribute to that end. It has been shown that even for a simplified case of one Taylor 
bubble passing through a heated tube, the hydrodynamics and their effect on heat transfer 
are quite complex. Consequently, one would expect that the model for this regime would 
need to be fairly advanced. It is pertinent to inspect previous models for their means of 
handling the complexities which may be incorporated into new, more thorough models. 
Several models have been developed to describe slug flow, including the 
aforementioned three-zone model by Thome et al. [40] for microchannels. The model 
breaks down an individual bubble’s unit cell into a liquid slug, the bubble, and a vapor slug 
which occurs if the liquid film around the bubble dries out. Given that the model was 
developed for capillary flow, the assumptions made of homogenous flow and uniform 
saturation conditions are reasonable for that situation. As has been shown, however, these 
assumptions are not valid for slug flow with Taylor bubbles which makes the extrapolation 


of the model to flow of this type less appropriate. 
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A model for macroscale, gas/liquid, vertical-upward slug flow was developed 
analytically by Barnea and Yacoub [35] in which an energy balance was performed on 
fluid cross-sections to determine the local wall and fluid temperatures as a function of time. 
They suggested that the transfer of heat from the wall to the fluid was largely dependent 
on the local velocity within the slug and film regions. The heat transfer coefficient in each 


was approximated using Colburn’s [68] correlation for forced convection inside tubes, 
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(12) 
where D,, is the hydraulic diameter (D in the slug, 26 in the film) and U is the mean liquid 
velocity for either the slug or film region. 

In an effort to predict peaks in heat transfer coefficient for near-zero vapor quality flows, 
Barbosa and Hewitt [69] developed a model for large scale turbulent liquid/vapor slug flow 
with Taylor bubbles. Their methodology utilized superposition to approximate the average 
heat transfer coefficients for the slug and film regions. In both cases, empirical correlations 
were superimposed to account for convective and nucleate boiling contributions. 

Within the slug, Barbosa and Hewitt suggested the Dittus-Boelter correlation for 
turbulent forced convection with a modified Reynolds number to account for entrainment 
of small bubbles. The calculated convective heat transfer coefficient was multiplied by 
Chen’s [54] two-phase enhancement factor, F. Heat transfer from nucleation was 
approximated with the Forster and Zuber [70] correlation multiplied by Chen’s [54] 
nucleate boiling suppression factor, S. In the absence of nucleation and entrained bubbles 


(as in the current experimental results), this formulation for the wake heat transfer provides 


a constant heat transfer coefficient dependent on the bulk liquid velocity within the slug. 
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The suggested prediction for the average convective film heat transfer was a correlation 
by Chun and Seban [71] for a turbulent falling film on the outside of a vertical heated 


cylinder, 
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Nucleation was again accounted for through Forster and Zuber’s correlation and Chen’s 
suppression factor. An interesting caveat relating to the film region was also included in 
the development of the model’s energy balance equations. The assumption was made that 
the liquid supplied to the film of a bubble by the leading liquid slug originates from an 
annular region located at the wall of the tube. A schematic of this theory (originally 
proposed by Orell and Rembrand [72]) is provided in Figure 3.23. The area at the center 
of the tube is forced upwards by the bubble, while the annular wall region flows downward 
forming the film. Barbosa and Hewitt extended this idea to the temperature field, 
hypothesizing that the average temperature within the film was higher than that of the slug 
because the cooler liquid from the core does not enter the film. 

It is evident from this summary of the currently available slug flow models that room 
for improvement exists in the prediction methods of heat transfer coefficient, especially in 
the wake where the mechanisms shown in this work are not addressed. While the models 
typically apply modified channel flow correlations, a more realistic approximation of 
plunging film heat transfer is the slot wall jet geometry (Figure 3.24) often seen in film 


cooling applications. 
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Figure 3.23: Schematic of slug liquid partitioning proposed by Orell and Rembrand [72]. 
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Figure 3.24: Schematic of slot injection geometry, adapted from Bittlinger et al. [73]. 


For the case of combustion chamber cooling, a thin film of cool air is injected via a slot 
near the chamber wall to protect it from the hot combustion gasses passing in the free 
stream. Several researchers have investigated the local heat transfer coefficient along a 


cooled wall in the presence of flowing heated air with various blowing ratios. The blowing 
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ratio, M, is defined as M = p¢U¢/PooUco, where pr is the film fluid density, Uy is the film 
velocity, Po. is the bulk fluid density, and U,, is the bulk fluid velocity. Bittlinger et al. [73] 
found that the heat transfer coefficient increased sharply just after the film injection and 
then decayed to a steady state value as the film diffused downstream as shown in Figure 
3.25. These profiles are remarkably similar to those found experimentally in Figure 3.16 
for rising Taylor bubbles. Bittlinger et al. noted that above a blowing ratio of 1.6 
(Us /U.>1) the peak in heat transfer coefficient increased with M. This finding is contrary 
to the current results, but it is hypothesized that the decrease in peak heat transfer 
coefficient for Taylor bubbles may be related to a transition of the background liquid flow 
from laminar to turbulent regime. The film cooling experiments were conducted under fully 
turbulent conditions. 

Seban [74] developed a correlation for local heat transfer coefficient downstream of the 


film injection for U;/U..>1 and a constant heat flux condition at the wall, 
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where s is the slot height. The heat transfer coefficient was defined as h = q"/(Ty — Taw), 
where Tyy is the adiabatic wall temperature (the temperature of the wall without heating 
or cooling given the same flow conditions). A comparison is made between Seban’s 
prediction for the modified heat transfer coefficient using experimental fluid properties and 
experimental results in Figure 3.26. While the prediction results are not directly 
comparable to experiments, it shows that a correlation of this type could be used to more 


realistically approximate the heat transfer profiles behind Taylor bubbles. Unfortunately, 
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due to the progression of film cooling research from slot-type to patterned hole geometries, 
recent predictive work on this topic is sparse. Future studies (whether experimental, 
numerical, or theoretical) on the flow characteristics and heat transfer of plunging jets may 


be required. 


s=3mm, t=imm 
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Figure 3.25: Heat transfer coefficient profiles for isothermal blowing T,, = Tr = 270°C for various 
blowing ratios from Bittlinger et al. [73] (used with permission). Here positive x indicates distance 
into the wake. 


The heat transfer prediction within the film region is also of importance, despite its 
relatively small effect on the overall heat removal for the cases described in this work. As 
mentioned in the review, modelers have proposed either modified single-phase internal 
flow correlations or falling film correlations for calculation of the heat transfer coefficient. 


A comparison of the suggested correlations to experimental data is shown in Figure 3.27 
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for a bubble rising in laminar background flow. It can be seen that Colburn’s correlation 
[68] for turbulent internal flow significantly over-predicts the heat transfer. Limitations of 
the experiment prevent measurement of the velocity profile within the film to determine if 
transition to turbulence occurs. In addition to Chun and Seban’s [71] correlation for a 
turbulent falling film, their correlation for a laminar falling film, 
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was also included in the figure for two cases. In the first case, the prediction was made 
using a constant film thickness and film velocity based on measurements of the fully 
developed film. This instance is representative of many model formulations where the 
bubble is assumed to be a cylindrical plug rather than bullet-shaped. The actual film 
thickness profile and calculated local film velocity were substituted in the second case. It 
is evident that the film heat transfer coefficient is also over-predicted by the turbulent 
falling film correlation, but agreement is seen between experimental results and the laminar 
film predictions. Furthermore, the inclusion of the local film thickness and velocity 
improved the agreement by better capturing the shape of the experimental curve. While it 
is unclear whether the mechanisms assumed within the laminar film correlation are directly 
applicable to the Taylor bubble film (especially if the theory of Barbosa and Hewitt 
regarding annular and core regions is considered), the prediction capability is evident and 
may serve useful in future model development. 

To conclude, it had been shown that particularly in the bubble wake, current slug flow 
models lack the ability to capture the magnitude and shape of the heat transfer coefficient 
profiles. Additional work will be required to apply slot film cooling mechanisms to the 
current results, but data from the field illustrates promise in that regard. Finally, it appears 
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that previously developed correlations for falling films may be directly applicable to 
modeling of the heat transfer coefficient for rising Taylor bubbles. 
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Figure 3.26: Comparison of predicted heat transfer coefficient by Seban [74] for slot injection film 
cooling to experimental results for Taylor bubble wake (experimental conditions: G=200 kg/m’s, 


a/g=1, q"=1700 W/m?, hsp=192 W/m’K). 
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Figure 3.27: Comparison of experimentally measured film heat transfer coefficient the correlations 
of Colburn [68] and Chun & Seban [71]. The conditions were: U; = 37 mm/s, U, = 142 mm/s, 


U, = 851 mm/s, Re, = 791, and q" = 1370 W/m’. 
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CHAPTER 4: BUBBLE EFFECT ON NUCLEATE BOILING 


The results presented in the previous chapter detailed the hydrodynamics and heat transfer 
associated with single elongated bubbles rising through a heated tube in which the initial 
condition was thermally-developing single-phase flow. While many applications exhibit 
these characteristics, some operating or emergency conditions in nuclear reactors result in 
nucleation at the wall. Bubbles nucleated near the bottom of the channel coalesce to form 
elongated bubbles and move upward, thereby affecting nucleation sites downstream. Using 
the same experimental setup described in Chapter 2, a preliminary study was completed on 
a single Taylor bubble rising through a tube where nucleation was present (see Figure 4.1). 
The conditions tested are summarized in Table 4.1. Power applied to the silicon tube was 
set such that nucleation was present over at least one half of the visible length to ensure 
data was collected over a sufficiently large area to be representative. Experiments were 
completed with three heat fluxes at each liquid mass flux to determine the influence of the 
bubble on increasingly vigorous nucleation. At elevated heat fluxes, the flow regime was 
observed to evolve significantly between the incipience of nucleation and the exit of the 
tube. Analyses for this study were isolated to the portion where nucleation was observed. 
The heat transfer within the silicon tube as the bubble passed was measured by the same 
technique described in Chapter 2. Due to the significant increase in heat transfer coefficient 
resulting from nucleation, thermal losses from the tube were approximated to be negligible 


at the end caps (conduction) and 5% to the ambient (natural convection). 
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Flow 


Figure 4.1: Schematic of initial conditions for single Taylor bubble influence on nucleate boiling 
tests. 


Table 4.1: Summary of nucleate boiling test parameters. 


G (kg/m’s) _ q" (W/m?) Tsub (CC) 


100 6630 2:/ 
9520 2.8 
14160 20 
185 9340 2.1 
13620 2 
18050 12 


As with the single-phase tests, the heat transfer coefficient was plotted as a function of 
distance from the bubble tail (z/D>0: ahead of tail, z/D<0: behind tail). Figure 4.2 shows 
the heat transfer coefficient profiles for G=100 kg/m*s and three heat fluxes. The 


measurements were left unscaled to emphasize the increase in overall heat transfer with 
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increasing heat flux. Unlike the results presented for the single-phase initial conditions in 
Chapter 3, the heat transfer coefficient was seen to decrease in the liquid film. This 
behavior can be explained if the ebullition model proposed by Hsu [75] is considered. He 
theorized that the activity of a wall nucleation site (growth of a vapor embryo) was limited 
by the growth of the thermal boundary layer, 6;, near the heated wall. If a sufficiently thick 
Or is not present, the thermal conditions around the embryo would not be conducive to 
nucleation. In the case of a passing Taylor bubble in nucleate flow boiling, the acceleration 
of the liquid near the wall due to film formation caused the thermal boundary layer to thin 
such that nucleation could no longer occur. This behavior was observed near the nose for 
all heat fluxes tested, with the magnitude of degradation increasing with q”. It is well 
known that the frequency of bubble nucleation and subsequently the heat transfer 
coefficient is related to the superheat of the heated surface. In this case, suppression of 
nucleation in the film became progressively more detrimental to the heat transfer 
coefficient as the heat flux, wall superheat, and agitation due to bubble departure increased. 
This can be visually confirmed from heat transfer contour plots similar to those presented 
in earlier sections. In this case, however, the motion of the bubble through the tube is 
observed from the laboratory reference frame as shown in Figure 4.3. A strip of low, 
relatively uniform heat transfer coefficient indicates the film region where heat transfer 
degradation occurred. Ahead and behind the strip/bubble, nucleation sites are visible as 
dots of high heat transfer. It is evident that an increase in applied heat flux increases the 


nucleation frequency and site density. 
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Figure 4.2: Heat transfer coefficient as a function of distance from the bubble tail as the bubble 
passes through nucleate boiling with G=100 kg/m?s and g"=6630, 9520, and 14160 W/m. 


The bubble length, and subsequently the distance ahead of the tail at which the heat 
transfer began to degrade, increased with heat flux due to higher evaporation rate of the 
bubbles as they passed through the silicon tube. As a consequence of the larger lengths and 
heat flux, re-nucleation within the liquid film was observed for the cases of q’=9520 and 
14160 W/m/?, leading to an increase in heat transfer coefficient prior to the bubble tail. 
Visual evidence of this behavior is provided in Figure 4.4 where a representative IR image 
illustrates small cold spots on the tube inner wall during the passage of the bubble which 


are indicative of nucleation sites. 
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Figure 4.3: Wake heat transfer contours for G=100 kg/m2s and heat fluxes: a) q"=6630 W/m2 and 
b) q”=9500 W/m2. 


The characteristics of suppression and re-nucleation observed in the liquid film are 
summarized in Figure 4.5 where process is divided into three zones (a, b, and c). In zone 
a, nucleate flow boiling is well established, resulting in steady heat transfer coefficient, 


wall temperatures, and thermal boundary layer thickness (6,). Formation of the film in zone 
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b accelerates the liquid and decreases 6; such that the superheat conditions are not sufficient 
for bubble growth. As a result, a decrease in the heat transfer coefficient is coupled with 
an increase in the inner wall temperature. Eventually, in zone c, the additional superheating 
of the wall overcomes the thinner boundary layer thickness to re-establish nucleation. The 
inner wall temperature was observed to stabilize slightly above (about 5%) the steady value 


ahead of the bubble due to the additional superheat requirement for nucleation. 
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Figure 4.4: Representative IR image showing the passing of a Taylor bubble where nucleation is 
present in the liquid film. The testing conditions were G=100 kg/m’s, q’=14160 W/m’, Tsup=2.1 
i OF 


The wake heat transfer enhancement due to vortex shedding (see Figure 4.2) was found 
to be significantly smaller than was observed with the single-phase initial conditions. This 
result is unsurprising given that the agitation caused by wake turbulence is likely 
comparable to forced nucleate boiling. The high levels of mixing also contribute to the 
rapid decay of the wake enhancement. Due to the absorption of the silicon and polyimide 


tape, details of the wake hydrodynamics could not be more closely examined. 
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Figure 4.5: Schematic of Taylor bubble rising through nucleate boiling with representative graphs 
of the measured heat transfer coefficient, wall temperatures, and approximated thermal boundary 
layer. 


Additional tests were completed at G=185 kg/m?s to determine the effect of mass flux 
on the heat transfer profiles. As can be seen in Figure 4.6, degradation of the heat transfer 
coefficient is still observed within the film due to suppression of nucleation. Re-nucleation 
is seen for both q”"=13620 and 18050 W/m”. It can also be pointed out that the heat transfer 
coefficient ahead of the bubbles at g"=9340 and 13620 W/m” is smaller than are seen in 


Figure 4.2 for similar heat fluxes of q”=9520 and 14160 W/m”. This can again be attributed 
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to a reduction of active nucleation sites caused by the thinner thermal boundary layer at a 


higher mass flux. 
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Figure 4.6: Heat transfer coefficient as a function of distance from the bubble tail as the bubble 
passes through nucleate boiling with G=185 kg/m?’s and q"=9340, 13620, and 18050 W/m7?. 


While preliminary, the results of this study provide insight into the effect of long 
coalesced bubbles on nucleate flow boiling. It is evident that the presence of sequential 
bubbles or bubble trains would be highly detrimental to a thermal management system 
operating under these conditions. Continued work on this topic can provide not only 
invaluable information to heat transfer modeling efforts, but also potential solutions to 


prevent elongated bubble formation and/or maintain the stability of nucleation sites. 
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CHAPTER 5: CONCLUSIONS 


The work presented in this dissertation is part of a broad effort by two-phase researchers 
to form accurate and reliable prediction methods for flow boiling heat transfer for earth- 
and space-based applications. To do so, previous research has indicated that two main 
topics are in need of further investigation: 1) development of flow boiling regime maps 
which account for fluid properties, heating conditions, surface properties, and gravity level; 
2) identification of prominent heat transfer mechanisms for each regime to be applied to 
heat transfer modeling. The focus of this study was to address the second area with an 
emphasis on the characteristics of slug flow. Specifically, individual unit cells of slug flow 
(single bubbles) were investigated so that the initial and boundary conditions could be more 


easily specified and the effect of the rising bubble on a heated tube identified. 


5.1 INTELLECTUAL CONTRIBUTIONS 


5.1.1 Adaptation of Flow Boiling Experiment 
The original UMD variable gravity flow boiling experiment was modified to increase its 
measurement capabilities and to create single elongated bubbles. The latter was 
accomplished by the addition of a bubble generation section consisting of two tubing 
segments connected downstream by a three way valve. One segment encased a wire heater 
used to vaporize liquid, creating a bubble of known volume. The second leg served as a 
bypass during bubble generation. 

An adiabatic section was also added downstream of the heat transfer measurement 
section to allow for measurement of the film thickness and observation of the bubble 


dynamics. Several film thickness sensors were tested before the Keyence LK-G5000 model 
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was chosen, purchased, and installed. A compact, high-speed Sentech camera was also 
purchased and installed to capture flow visualization videos. 

Finally, the experiment was prepared for variable gravity experiments aboard NASA’s 
DC-9 aircraft. This process included verification of the experiment’s structural 
requirements, completion of NASA safety documentation, and transport of the rig to 
Ellington Field in Texas. Successful flight testing was completed, the results of which have 


been presented in the previous chapters. 


5.1.2. Characterization of Elongated Bubbles 

A study on the hydrodynamics and heat transfer in the presence of an elongated bubble 
rising co-currently with vertical, upward flowing liquid was completed. Data collected in 
the laboratory and on an aircraft flying parabolic trajectories allowed for bubbles 
characterized by a wide range of Bond numbers to be observed. Visual analysis of the 
bubbles illustrated the effect of gravity, and subsequently drift velocity on the bubble 
profile. Bubbles exhibited small drift velocities and shapes characteristic of capillary flows 
at low gravity levels, while terrestrial and hypergravity conditions produced Taylor bubbles 
moving at higher drift velocities. The measured film thicknesses were shown to vary little 
with acceleration and be in agreement with correlations. 

The flow dynamics induced by the different gravity levels was found to strongly 
influence the bubble shape and the local heat transfer characteristics. Large heat transfer 
enhancement was observed behind Taylor bubbles due to the presence of vortices in the 
wake. Capillary bubbles provided little wake heat transfer enhancement and only slight 
improvement in the liquid film. The frequency of vortex shedding for Taylor bubbles was 


shown to increase linearly with the plunging velocity, allowing for the proposal of a critical 
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impinging velocity to predict vortex incipience. Wall jet separation was found to affect the 
tube wall at differing penetration lengths depending on the acceleration and bubble drift 
velocity. Time-averaged heat transfer coefficient profiles showed that the variation in 
penetration length coincides with smaller and broader peaks in wake enhancement for 
a/g=1.8 compared to a/g=1. 

Comparisons made with numerical simulations completed by Dr. Mirco Magnini at 
EPFL, Switzerland showed strong agreement in the bubble hydrodynamics as well as heat 
transfer characteristics. The local flow and temperature fields provided by the simulations 
enabled a better understanding of the vortex shedding process in the wake and its effect on 
the heat transfer coefficient. 

The information gleaned from experiments and simulations were then applied to slug 
flow heat transfer models in the literature to determine their validity. It was found that the 
wake heat transfer coefficient was typically approximated using turbulent single-phase 
correlations. A more appropriate representation of the wake dynamics was proposed to be 
the slot film geometry, which produces similar heat transfer profiles to Taylor bubbles. 
Heat transfer within the film was shown to be reasonably predicted by laminar falling film 


correlations. 


5.1.3 Taylor Bubble Effect on Nucleate Flow Boiling 

A preliminary study was completed on the effect of a single Taylor bubble rising through 
a heated tube in which nucleate boiling had been established. Results showed that the 
acceleration of liquid into the bubble film thinned the thermal boundary layer and 
suppressed nucleation, leading to a decrease in heat transfer coefficient. The effect of 


suppression was more pronounced (larger degradation) at higher heat fluxes as more 
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nucleation sites were active prior to the bubble’s arrival. In cases of long bubbles and high 
heat flux, the thermal boundary layer grew quickly within the film allowing for nucleation 


to be re-established. 


5.2 RECOMMENDATIONS FOR FUTURE WORK 


The results of this study, while providing valuable insight into the behavior and heat 
transfer around elongated bubbles, illustrate the complexity of slug flow and the need for 
further research in this area. Several directions could be taken to improve and advance 


upon the current work: 


1. Given the insight provided by the presented experimental results and the overall goal 
of this study, a logical step in the progression of this work would be to develop a model 
for a single Taylor bubble rising in co-current flow through a heated tube. General 
considerations for heat transfer in the wake region were discussed in Section 3.4, where 
treatment of the plunging jet as a wall jet was suggested. It was also shown that the film 
heat transfer could be reasonably approximated using correlations for a falling film. A 
successful model for single Taylor bubbles can then be expanded to bubble trains, 
which are of more practical interest. 

2. Improvements to the current experiment could provide a more detailed and complete 
data set allowing for the resolution of several outstanding questions. Replacement of 
the silicon tube with a sapphire tube of the same size would allow for simultaneous 
heat transfer measurements, flow visualization, and film thickness measurements to be 
made within the heated section. This would improve greatly the characterization of the 


tail motion effect on wake heat transfer as well as the variation in film thickness due to 
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evaporation. An additional, and more complicated, modification would be the addition 
of particle image velocitmetry (PIV) or particle shadow velocimetry (PSV). Tracking 
of particles in the bubble wake would provide local, time-resolved visualization of the 
wake flow patterns dictating heat transfer at the wall. 

. As the flow field and heat transfer around single elongated bubbles has been relatively 
well studied, the current (and potentially improved) experiment could be applied to 
investigate the effect of sequential bubbles and bubble trains on the local heat transfer. 
Some work on sequential bubbles is currently being completed by the group of Shemer 
and Barnea at Tel-Aviv University [76], but the combination of local heat transfer 
measurements and flow visualization is lacking. A numerical investigation by Magnini 
et al. [42] showed that for two bubbles flowing through a microchannel, the second 
bubble caused larger heat transfer enhancement than the first. This was attributed to the 
overlap of the hydrodynamic effects on wall temperature. It is expected that similar 
findings would be obtained for Taylor bubbles in larger channels. 

. The preliminary study on Taylor bubbles rising through a heated tube where nucleate 
boiling is present illustrated that, unlike previous tests with single-phase flow, the 
bubble hydrodynamics were detrimental to heat transfer due to suppression of 
nucleation. Using an improved sapphire test section where nucleation sites can be more 
easily observed, the fundamentals of this behavior can be closely investigated. 

. Amore practical line of work may also be grown out of the early nucleation findings. 
If this type of flow (nucleation with passing large bubbles) were to occur in an 
application where heat transfer maximization were desired, it would be advantageous 


to either break up the large bubbles or to ensure that nucleation sites remained active 
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even in the presence of the bubble. A solution to the latter could be modification of the 
surface roughness of the tube wall such that a thinner boundary layer would be required 
for nucleation. This may have the additional effect of promoting the incipience of 
boiling at a lower superheat, which is advantageous from a heat transfer perspective. A 
study on surface modification of flow boiling channels could investigate these 


propositions. 
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APPENDIX A: CALIBRATIONS AND UNCERTAINTY 


A.1 UNCERTAINTY OVERVIEW 


The uncertainty analyses of this experiment were based on the practices laid out in 
Beckwith et al. [77] and Taylor [56]. The overall uncertainty of a measurement can be 
broken into two components: the bias uncertainty and random uncertainty, as shown in 
Figure A.1. Bias uncertainty defines the possible measurement error as an offset with 
respect to the actual value, while random uncertainty describes possible error in the 
measurement due to randomness in a signal or repeated measurements. Each of these 


uncertainty types can be calculated separately and then combined, as will be shown below. 


—— Signal 
Actual Value 
Mean Value 


Bias Uncertainty 
Random Uncertainty 


Dependent Variable 


Independent Variable 


Figure A.1: Description of measurement uncertainty categories. 
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A.1.1 Bias Uncertainty 

The bias uncertainty for measurements obtained during this study, particularly heat flux 
and heat transfer coefficients, were determined on a case-by-case basis as the calculation 
procedures vary. These analyses will be described individually in the sections that follow. 
One common method in all cases is accounting for propagation of error due to secondary 
or tertiary measurements used in the primary measurement calculation. For example, a 
measured value y is dependent on properties x1, *2,...,Xp,. The uncertainty of y, dy, is 


calculated by the equation, 


by = (Zox,) + (25x) foe (2ox,) (16) 


x4 
where, 6X1, 0X2, ..., OX, are the uncertainties of properties x1, Xz, ..., Xp». This equation may 
also be applied when combining the bias and random uncertainties to obtain the overall 


uncertainty of a measurement. 


A.1.2. Random Uncertainty 
Random uncertainty of calculations, calibrations, or transducer measurements was 
determined using the Student’s t-distribution method, where the distribution of a quantity 


t is defined by, 


pa 
Sx/VN 


(17) 


Here, x is the mean value for a sample, y is the mean value for a population, S, is the 
standard deviation of a sample, and N is the number of samples. The values for t are widely 
available in tabulated form and are dependent on the number of degrees of freedom (v = 
n—1) and the desired confidence, c. For all uncertainty calculations presented, the 


confidence was chosen to be 90%. The random uncertainty is defined by the quantity x — 
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tS). 


h= which represents the deviation of the sample mean from the population mean 


caused by relatively few measurements. 


A.1.3 Curve Fitting Uncertainty 

A final category of uncertainty calculation is the potential error induced by fitting a linear 
curve to acquired data, as was the case for calibrations described below. A method 
presented by Taylor [56] was utilized to predict the induced slope (B) and offset (A) 


uncertainty, 


ss yx? 

Oa = Cy eet oe 
N 

99 = Oy laa OaF 


where x is the independent variable, y is the dependent variable, and N is the number of 
samples. The parameter oy is an error estimation of the collected data with respect to the 


fitted line and is defined as, 


le 
N-2 


Oy = N 4 - A Bx)? (20) 


A.2 INDIVIDUAL MEASUREMENT ANALYSES 


A.2.1 Temperature 

Proper measurement of temperature is essential in heat transfer experiments to reduce the 
uncertainty of results and to ensure that the desired testing conditions are obtained. This 
experiment utilizes T-type thermocouples and an Electrophysics 660M infrared camera to 
measure the various temperatures in the system. Thermocouples were calibrated by 


immersion in a controlled water bath, with the reference temperature measured using a 
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NIST traceable mercury thermometer accurate to 0.1°C. The most important thermocouple 
measurements for experiments were the temperature at the test section inlet (Tin) and the 
temperature at the preheater inlet (Tpx,in); for optical property measurements and camera 
calibration the most important were the temperature of the blackbody (Th») and the ambient 
temperature within the containment chamber (Topbox). During thermocouple calibration, nine 
measurements were taken and compared to the reference thermometer as shown in Figure 
A.2a. Linear fits were computed for each thermocouple to obtain a calibration curve, which 
was used in data processing. The uncertainties of the thermocouples were calculated and 
include the error associated with the linear fit as well as the uncertainty of the mercury 


thermometer. A summary of the uncertainties can be found in Table A.1. 
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Figure A.2: Calibration curves for a) experiment and calibration thermocouples, b) infrared camera. 


A similar calibration procedure was conducted for the infrared camera. A blackbody 
(described in Kim et al. [47]) was heated and the opening observed by the camera. The 


measured camera temperatures were plotted against the measured blackbody temperature 
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as shown in Figure A.2b. The uncertainty of the camera measurement was found to be 


0.14°C and included the uncertainties of the calibration curve and blackbody temperature. 


Table A.1: Summary of important temperature uncertainties. 


Temperature Uncertainty (°C) Typical Value (°C) 


ie 0.13 58 
Tatas: 0.12 34 
Ts 0.12 30-90 
Tes 0.39 25 
Ties 0.14 50-70 


A.2.2 Flow 

The flowmeter used in this experiment, an Omega FLR1009, was calibrated by flowing 
water at constant rates through the flowmeter and collecting samples in a graduated 
cylinder over a specified duration. The graduated cylinder had a capacity of 100 mL and 
graduations of 1 mL. A stopwatch with a readout to 0.1s was used to record the time. Figure 
A.3 illustrates the measured flow rates as a function of output voltage. A linear curve was 
fit to the data to create a calibration curve used in the LabVIEW™ VI and post-processing. 
The uncertainty in determining the flow rate with the calibration curve was found to be 1.3 


mL/min (typical value 50-250 mL/min). 
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Figure A.3: Calibration data and curve for FLR1009 flowmeter. 


Of greater use in experimental analysis, however, is the liquid mass flux, rather than the 
flow rate. The mass flux can be determined from the flow rate with knowledge of the liquid 
density (p,) and tube area (A;) by equation (21), where the flow rate (Q) is now in the units 


of m?/s. 


— Qpi 
G=7- (21) 


The uncertainty of the mass flux is therefore dependent on the uncertainty of the flow 
rate, the liquid density, and the radius of the tube. As was discussed above, the flow rate is 
known to within 1.3 mL/min (2.23x10° m?/s). The uncertainty in HFE 7100 density will 
be determined in later sections, but was found to be 0.45 kg/m?. Finally, the tube radius 
was measured using precision calipers to within 0.01 mm. Using the propagation of 
uncertainties method, the overall mass flux uncertainty was found to be 3.95 kg/m°s 


(typical value of 50-200 kg/m’s). 
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A.2.3 Pressure 

The absolute pressure at the inlet of the heated silicon test section was measured using an 
Omega PX209-030A5V transducer with a range of 0-30 psia. The factory calibration 
(NIST traceable) was used as a baseline, but checks were completed using a custom 
manometer built in the laboratory. Figure A.4 shows the applied absolute pressure as a 
function of transducer output voltage. The calculated uncertainty of the calibration curve 


was found to be 0.01 psi. 
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Figure A.4: Absolute pressure transducer calibration. 
A.2.4 Fluid Properties 
Properties for liquid and vapor HFE 7100 were determined using curve fits to data provided 
by 3M (the manufacturer of HFE 7100). The absolute pressure at the test section inlet was 
used as the independent variable for property determination during testing and data analysis 
of saturated conditions. Local temperature, measured with T-type thermocouples, was used 
as the independent variable for single-phase property calculations. Uncertainties of the 


fluid properties were calculated by accounting for the error associated with the curve fit to 
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the data provided by 3M as well as the uncertainty in measurement of the independent 


variable. A summary of the major uncertainties is given in Table A.2. 


Table A.2: Summary of fluid property uncertainties. 


Property Uncertainty Typical Value 
Tsat PC] 0.14 60 

p, [kg/m?] 0.45 1390 

hy [MJ/kg] 0.97 112 

ut, [cP] 0.001 0.375 


A.2.5 Film Thickness 

The film thickness sensor (Keyence LK-G5000) utilizes a triangulation method of 
measurement as illustrated in Figure A.5. Laser light exits the sensor head at a specified 
angle and is reflected from the inner and outer walls of the tube and from the bubble 
liquid/vapor interface. Depending on the indices of refraction and thicknesses of the 
various layers, the reflected beams will intercept the sensor measurement array at different 
locations. The location at which the light reflects from the liquid/vapor interface is 


dependent on the film thickness as well as the inclination of the interface. 
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Figure A.5: Schematic of triangulation technique used by film thickness sensor. 


Assuming that the liquid film is smooth and the film thickness is not changing (the 
liquid/vapor interface is parallel to the tube wall), the multilayer can be approximated using 


flat plates to simplify the setup for calibration (Figure A.6). 


HFE 7100 


Microscope Glass 


Polyimide 


Fused Silica 


Figure A.6: Schematic of parallel calibration setup. 


A fused silica wafer with a thickness of 1.03 mm and index of refraction of 1.456 was 
used to simulate the Pyrex tube of 1.00 mm thickness and index of refraction of 1.474. 
Polyimide strips of known thickness (48+1 jam and 24+1 jim) were used as spacers between 


the silica wafer and a microscope cover class that simulated the bubble interface. The gap 
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created by the spacers was filled with HFE 7100 to replicate the bubble liquid film. The 


known film thickness as a function of the sensor readout was obtained over the range of 


film thicknesses measured in our tests and a linear fit to the data was obtained as shown in 


Figure A.7a. Under these conditions, the thickness could be measured to a minimum 


accuracy of 16 jm. 
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Figure A.7: Calibration curve for Keyence film thickness sensor: a) flat configuration and b) angled 
configuration. 


Shorter bubbles possessed developing liquid films which were not parallel to the tube 


wall. Characterizing the uncertainty associated with this angled film was completed in a 


similar manner to the original calibration. The polyimide films were unevenly stacked such 


that the microscope glass was inclined relative to the silica wafer as shown in Figure A.8. 
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Figure A.8: Schematic of angled calibration setup. 


Measurements were obtained at various inclination angles up to 0.27° (maximum 
bubble film angle measured was 0.26°) and the error relative to the actual film height was 
calculated. The error in measured thickness increased with increasing film angle as shown 
in Figure A.7b as expected. The maximum uncertainty for film thickness with the added 
correction required for the non-parallel liquid/vapor interfaces was 17 ym. 

During parabolic flight testing, settling of the experiment due to gravity variations 
yielded the Keyence sensor inoperable. Film thickness measurements were instead made 
using detailed image analysis of the high-speed video. In order to correct for optical 
distortion caused by the curvature of the glass tube, a calibration was completed in which 
a grid with 0.5 mm spacing was inserted into the tube and submerged in liquid HFE 7100 
as shown in Figure A.9a. The radial distance from the tube centerline from which each grid 
line was observed was correlated to the known grid spacing to create a calibration curve 
shown in Figure A.9b. Estimating the uncertainty for this technique includes several 
components and assumptions. First, the grid spacing is assumed to be accurate (it was 
printed using a 600 dpi laser printer and was checked using a ruler before insertion into the 


tube). Second, the potential error in determining the center each grid line was assumed to 
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be +1 pixel (+37 pm), yielding an uncertainty of = = 19 ym. This value was included 


with the calculated uncertainty of linear curve fitting to obtain a total uncertainty of 25 pm. 
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Figure A.9: Calibration of visual film thickness measurements: a) 0.5 mm grid inserted in tube and 
submerged with HFE 7100, b) calibration curve relating raw radial position to known grid 
distances. 


Film thickness measurements using the two methods were compared to ensure 
consistency and provide a means of validation. Bubbles rising in liquid with drift velocities 
over the range 105 < Ug < 163 mm/s were measured by both techniques. The film 
thickness at the bubble tail was plotted as a function of bubble length in Figure A.10. 
Measurement types are seen to agree with each other as well as with correlations by Brown 


[19] and Llewellin et al. [26] for fully developed films. 
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Figure A.10: Comparison of film thicknesses measured by a Keyence laser displacement sensor 
and analysis of high-speed visual images. 


A.2.6 Optical Properties 
The infrared technique utilized in this study for heat transfer measurements was dependent 
on knowledge of the optical properties of the materials constructing the tube multilayer: 
the silicon tube, Polyimide tape, and black paint. A rigorous series of experiments was 
conducted to obtain the absorptivity, reflectivity, and emissivity of the multilayer 
components according to the procedure described in Kim et al. [47], which is summarized 
below. 

The silicon tube wall and Polyimide tape were considered to be plane walls for the 
infrared technique’s optical model as the width of the camera’s pixel (0.11 mm) is 


significantly smaller than the diameter of the tube (6 mm). Using this assumption, ray- 
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tracing of radiation incident on a single plane wall yields equations for the apparent 


reflectivity (eqn. (22)) and apparent transmissivity (eqn. (23)) of the wall, 


= Day Pata 


ee a 7 eee 


2 
-_ (1 a Pree) Tm 
Tapp,m—o ~ | 2 Picts? ie 
- m 


(22) 


(23) 


where T,, and Pm—o are the transmissivity of the reflectivity of the wall, respectively. In 


order to extract the two absolute properties from these equations, two experiments are 


required to determine the apparent properties and then the equations may be solved 


simultaneously. The experimental setups are shown below in Figure A.11. 
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Figure A.11: Experimental setup for multi-layer optical property determination: a) apparent 


transmissivity measurement, b) apparent reflectivity measurement. 
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The energy received by the camera is a sum of the various sources, namely, the 
blackbody, emission of the film, and emission of the ambient surroundings. Equations for 
the energy balances are given below for the apparent transmission test with/without the 
blackbody (Figure A.11a), and for the apparent reflection test with/without the blackbody 
(Figure A.11b), respectively: 


Transmission: 


= 4 
Ez1 = Tapp,meo(Fi,-2,9Tp*) 


(24) 
+ Eapp,m—oo cones OT m) + Papp,m—oo (Piya, OT." ) 
Ex2 = Tapp,m—co (Fi,-a, OT co ) 
(25) 
+ Eapp,m—oo (ie OT mn) + Papp,m—oo (Piya, OT.) 
Reflection: 
Ey = Tapp,m—co (Fi,-a, OT x") 
(26) 
+ Eapp m—c(Fa,-a29Tm) + Pappm-«(Fi,-229Tb) 
E.2 = Tapp,mo — (Fi,-a, ee) 
(27) 


+ Eapp,m-oo (Cre OT mn") + Papp,m-oo (Pens oT") 


where Egppm—co is the apparent emission of the film/substrate and Fy,-,, is the fractional 
function defining the portion of Planck’s function absorbed by the camera. Subtraction of 
eqn. (25) from egn. (24) and eqn. (27) from egn. (26) yields simplified relations for the 
apparent transmission and reflection of the film/substrate, which were used to 


simultaneously solve eqns. (22) and (23) for the actual transmissivity and reflectivity. 
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Tapp,m-o = 
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= (29) 
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Papp,m-o = 


Uncertainties of the transmissivity and reflectivity of the polyimide tape and silicon tube 
were estimated by performing a perturbation analysis on the property determination 
procedure. This method entailed varying each of the inputs to determine their relative 
impact on the overall uncertainty. The major contributors were found to be temperatures 
measured by the camera and thermocouples, with uncertainties of 0.14°C and 0.12°C, 
respectively. Combining the influences of each input yielded the uncertainties summarized 
in Table A.3. The transmissivity of both silicon and polyimide tape was converted to the 


_ In(tm) 


form of absorption coefficient (defined as a = , where d is the sample thickness) 


to be applied in heat transfer calculations. By doing so, the uncertainty in thickness of the 
samples tested was introduced. Scanning electron microscope (SEM) images were used to 
obtain high-resolution thickness measurements of the polyimide tape, while a dial 
micrometer was used for the silicon tube. Uncertainties in thickness of the polyimide tape 
and silicon tube were found to be 0.32 jim and 10 pm, respectively. Values of absorption 


coefficient are also included in Table A.3. 
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Table A.3: Uncertainties associated with optical property determination. 


Parameter Uncertainty Typical Value 
Tsi[-] 0.0064 0.9 

TpI-] 0.0072 0.74 

psi l-] 0.0026 0.3 

Pp I-] 0.0019 0.044 

asi [m*] 6.5 103 

ap [m*] 167 B27) 


A.2.7 Effective T.. Calibration and Uncertainty 

As described in Chapter 2, an in-situ calibration was completed during testing to determine 
an effective T.. to be used in the heat flux calculations. This calibration removed the effect 
of reflections within the tube, variation in optical properties due to thermal gradients within 
the multilayer, and non-uniformity in the surrounding temperature due to heating of the 
camera housing. The tube wall temperatures were measured as single-phase liquid flowed 
through the silicon tube, which was heated with a constant power input. In order to calibrate 
the effective surroundings temperature, knowledge of the heat transferred into the fluid was 
required. Losses from the tube included natural convection to the surroundings and 
conduction into the connections at each end, each of which was estimated using the outer 
tube wall temperature profile (an example is shown in Figure A.12). Conduction losses 
were approximated by measuring the thermal gradient at the outlet of the silicon tube. The 
temperature gradient at the inlet was not used for these calculations as the profile was 
dominated by thermal entry region effects. Instead, the inlet losses were assumed to be the 
same as the outlet given that the O-ring connections were identical. Including both ends, 


conduction losses were approximated to be Qiosscaps = 23+3%. 
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Losses to the surroundings by natural convection were estimated using a correlation by 
Churchill and Chu [49] for free convection from a vertical wall, 

Nu, = 0.68 + 0.515Ra,/* (30) 
where Nuy is the wall-averaged Nusselt number and Ray is the wall averaged Rayleigh 
number (based on the average wall temperature - calculated from the profile). Loss to the 
surroundings was approximated to be Qigss conv = 941%. Together, the total power loss 
from the tube was Qiossr = 324%, which compare well with the fluid energy balance 


calculation of Qigss.rc = 35412%. 
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Figure A.12: Representative wall temperature profile showing thermal gradients at tube 
connections. 


The net power entering the fluid was converted to heat flux using the geometric 
parameters of the tube and applied as an input to calculate the effective T.. from the data 
reduction MATLAB scripts. The uncertainty in the surroundings temperature was 


estimated to be +0.2°C (+1%), accounting for potential error in the heat loss calculations. 
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A.2.8 Heat Transfer 

The heat flux and heat transfer coefficient procedure outlined in Chapter 2 was carried out 
within a MATLAB script due to the complexity and breadth of the calculations. 
Uncertainty of these measurements was approximated using a perturbation analysis in 
which inputs to the calculations were perturbed slightly from the assumed value to 
determine the impact of each input’s uncertainty on the overall uncertainty. The heat flux 
calculation was a function of several parameters related to the IR temperature measurement 


and the coupled optical/inverse conduction problem, 
q" = f Tei» Teo, ksi, Kp, Kaw dsj, dy, daa, Psi-co» Psi-p» sis Cp, Ep, To) 


where T,; and T,, are the inner and outer wall temperatures measured by the IR camera; 
Ks, kp, and Kgq are the thermal conductivities of the silicon, polyimide, and adhesive; ds;, 
d,,, and dgq are the thicknesses of the silicon, polyimide, and adhesive layers; psj_.. and 
Psi-p are the reflectivities of silicon to air and silicon to polyimide; a; and @, are the 
absorption coefficients for silicon and polyimide; €, is the emissivity of the black paint; 
T.. is the effective temperature of the surroundings based on calibrations. A representative 


data set was selected and each input was independently varied by +10% and +20% to create 


" 


a curve from which the slope ($*) could be calculated at the assumed value. The slopes 


were then used to calculate the overall bias uncertainty via the propagation of error method 
(see eqn. (16)) with the corresponding input uncertainties. Major contributors to the bias 
uncertainty were found to be T;;, Psi, Psi-p» %si, and @,. This procedure was repeated 


for six values of measured heat flux to determine if the uncertainty is a function of the 


116 


measurement itself. Figure A.13 illustrates that the bias heat flux uncertainty varies little 
over the range of heat fluxes tested, with percentages ranging from 110% at q" = 0.2 


kW/m? down to 20% at q" = 11.7 kW/m”. 
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Figure A.13: Bias uncertainty in heat flux as a function of measured heat flux. 


Estimation of the bias uncertainty in heat transfer coefficient was completed using a 
procedure similar to that used for the heat flux. Specifically, the local heat transfer 


coefficient was defined as 


n= Gta) (31) 


where q" and T,; are calculated from the MATLAB script. As before, a perturbation 
analysis was completed to determine the relative impact of the various input parameters at 
several values of the heat transfer coefficient. Final values of heat transfer coefficient bias 


uncertainty are plotted in Figure A.14 as a function of heat transfer coefficient and include 
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contributions of q", T;;, and T,,. The bias uncertainty was found to increase significantly 
with increasing heat transfer coefficient from 302 at h=266 W/mK (113%) to 733 at 
h=2315 W/m°K (32%). These large uncertainties were partially negated in the presentation 


of results by offsetting h with hsp, as described in the IR technique description. 
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Figure A.14: Bias uncertainty in heat transfer coefficient as a function of measured heat transfer 
coefficient. 
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APPENDIX B: FILM VELOCITY DERIVATION 


In order to determine the velocity of the liquid film moving toward the rear of the bubble, 
we employ a conservation of volumetric flow rate. This is done by setting the flow rate 
calculated at a control surface in front of the bubble equal to the flow rate at a control 


surface in the film. 


R R 
C= | [U,(r) — Up] + 2ardr = | [U, + Up] - 2ardr 
0 R-6 


R R 
RU, — | U(r) + 2nrdr = -| [U, + Up]: 2ardr 
0 R-6 


If we want the average film velocity at each z location then U; is a function of z only. 


Additionally, U;, is a constant and can either be predicted or experimentally measured. The 


right hand side is given by, 
R 
— | [Up + Up] 2urdr = —n(U; + Up)(2R6 — 5?) 
R-6 


Combine with the left hand side and simplify to obtain, 
2 i U,(r)rdr — U,R? 
a (2R65 — 52) = 
By integrating over the liquid velocity profile, equation (1) can be obtained. 


(U, — Ug) R* 


Ui = ORs = 62 
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